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by 

* 
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This  report  contains  programs  used  in  computing  the  tables  presented 
in  Jerome  Bracken  and  Arthur  Schleifer,  Jr.,  Tables  for  Normal  Sampling  with 
Unknown  Variance;  The  Student  Distribution  and  Economically  Optimal  Sampling 
Plans ,  Division  of  Research,  Graduate  School  of  Business  Administration,  Har¬ 
vard  University,  196^.  The  tables  are  essentially  of  two  kinds:  tables  of 
the  ordinary  Student  ntM  density  and  cumulative  functions,  and  tables  to  facil¬ 
itate  Bayesian  analysis  of  certain  commonly  occurring  decision  problems  in 
which  sampling  may  or  may  not  be  involved.  The  programs  given  in  this  report 
could  be  used  to  compute  tables  for  parameter  values  other  than  those  of  the 
book  either  by  reading  in  alternative  data,  or  by  straightforward  modification 
where  the  parameters  of  the  book  are  included  in  the  programs.  The  programs 
are  written  in  FORTRAN  II,  and  in  computing  the  tables  in  the  book  they  were 
used  on  the  IBM  7090  and  IBM  1401.  It  should  be  noted  that  new  programs  have 
been  written  to  perform  some  of  the  computations  faster,  more  accurately,  or 
more  efficiently. 
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1 .  Introduction 


This  report  contains  the  computer  programs  used  in  computing  the  tables 
presented  in  Bracken  and  Schleifer  [l].  The  tables  are  essentially  of  two  kinds 
tables  of  Student  "tn  density  and  cumulative  functions,  and  tables  to  facilitate 
Bayesian  analysis  of  certain  commonly  occurring  decision  problems  in  which  sam¬ 
pling  may  or  may  not  be  involved.  Theoretical  aspects  of  the  Bayesian  analysis 
and  proofs  of  the  formulas  used  in  calculating  the  Bayesian  tables  are  given  in 
Schleifer  [3 ],  and  extend  results  given  in  Raiffa  and  Schlaifer  [2]. 

Bracken  and  Schleifer  [l ]  includes  an  introduction  which  describes  the 
tables,  discusses  prior,  posterior  and  preposterior  properties  of  Student  dis¬ 
tributions,  presents  the  two-action  problems  with  linear  utility  which  can  be 
analyzed  by  means  of  the  optimal  sample  size  and  maximal  net  gain  tables,  and 
discusses  other  applications  of  the  tables.  Interpolation  in  the  tables  and 
details  of  their  computation  and  verification  are  considered. 

Although  other  tables  of  the  Student  "t"  density  and  cumulative  functions 
have  been  published,  the  most  notable  being  Smirnov  [4],  the  computer  programs 
used  in  their  computation  have  not  been  published  or  described. 
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2.  Computational  Methods  and  Computers  Used 

Section  6  of  Bracken  and  Schleifer  [l  ]  discusses  computation  and  verifica¬ 
tion  of  the  tables  in  some  detail,  and  should  be  used  in  conjunction  with  the 
programs  given  in  this  report.  The  programs  are  organized  in  a  manner  which 
closely  parallels  the  discussion  of  their  computation  and  verification. 

All  of  the  computer  programs  are  in  FORTRAN  II,  and  were  written  for  one 
or  the  other  of  two  IBM  computers  --  the  709*9  anc*  the  1401.  The  tables  of  the 
Student  density,  right-tail,  and  linear-loss  functions  fg^  ,  Gg^  and 
were  computed  on  the  7090  in  double  precision.  They  were  further  edited  for 
printing  on  the  7090  in  single  precision.  Included  as  starting  information  in 
the  double  precision  computations  were  values  of  B(^,^v)  and  Gg^.(t=10|v)  to 
1 6  significant  figures,  which  had  previously  been  computed  on  the  1401.  It 
should  be  noted  that  a  program  for  computing  has  been  written  in  SPS 

for  the  IBM  1620,  but  is  not  included  here.  The  optimal  sample  size  and  net 
gain  functions  were  computed  on  the  7090  in  double  precision.  Supplementary 
programs  for  computing  optimal  sample  size  and  net  gain,  written  in  FORTRAN  II 
for  the  1401  and  used  to  check  the  computations  on  the  7090,  are  also  included 
in  this  report. 

3 .  Computer  Program  Summary 

This  section  presents  a  summary  of  the  computer  programs.  Inputs  and 
outputs  of  each  program  are  given.  Programs  of  sections  A,  B,  D,  and  E  were 
run  on  the  7090,  programs  of  sections  C  and  F  on  the  1401. 

A.  Programs  for  Computing  Tables  fg^  ,  and  Lg^ 

A.1  Standardized  Student  Functions  for  v  <  00 

Inputs  -  none  (contains  Gg^(t=10|v)  and  B(-§,Jrv)  for  all  v's) 
Outputs  -  fg*  ,  Gg*  and  ,  t  =  0(.0l)l0 

A. 2  Standardized  Student  Functions  for  v  =  »  (Normal) 

Inputs  -  none  (contains  Gg^.( t=1 0 1  v=oo) ) 

Outputs  -  fg^.  ,  Gg^.  and  Lg^.  ,  t  =  0(.0l)l0 
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B.  Program  for  Printing  Tables  ,  Gg^  and 

Inputs  -  fg^  ,  and  Lg^  ,  t  =  0(.0l)l0  ,  double  line  indicator 

Outputs  -  fg^  ,  Gg^.  and  Lg^  tables  in  final  form 

C.  Supplementary  Programs  for  Computing  Standardized  Student  Functions 
C.1  Evaluation  of 

Inputs  -  v  and  values  of  Riemann!s  zeta  function  for  integral  arguments 
Outputs  - 

C.2  Finite  Series  Expansion  for  Gg*.(t|v)  >  Integral  Even  v  <  00 
Inputs  -  v  ,  B(f,-g-v)  y  t 

Outputs  -  ^g*(tlv)  >  Gg*(tlv)  >  Ls*(tlv) 

C-3  Finite  Series  Expansion  for  Gg^(t|v)  P  Integral  Odd  v  <  00 
Inputs  -  v  y  B(|^-|v)  ,  t 
Outputs  -  fg*(t|v)  ,  Gg*(t|v)  ,  Lg^(t|v) 

i  2 

C.4  Infinite  Series  Expansion  for  Gg*(t|v)  }  Real  v  <  00  }  t  <  v 

Inputs  -  v  y  B(|,|v)  ,  t 
Outputs  -  fg*(t|v)  ,  Gg#(t|v)  ,  Lg*(t|v) 

i  2 

C.5  Infinite  Series  Expansion  for  Gg^(t|v)  y  Real  v  <  00  y  t  >  v 

Inputs  -  v  y  B(|-,tv)  ,  t 
Outputs  -  fg*(t|v)  ,  Gg#(t|v)  ,  Lg*(t | v) 

C.6  Continued  Fraction  for  Gg^(t|v=°°)  ,  t  >  2 
Inputs  -  t  y  n  (number  of  terms  in  series) 

Outputs  -  Gg^(t|v=°°) 

C. 7  Series  Expansion  for  Gg^(t|v=°°)  ,  t  <  2 

Inputs  -  t  ,  n  (maximum  number  of  terms  in  series) 

Outputs  -  Gg^.(t|v=00) 

D.  Programs  for  Computing  and  Printing  Optimal  Sample  Size  (R)  and  Maximal 
Net  Gain  (r) 

D.  1  R  and  T  for  v  «  1^,2 

Inputs  -  index  on  v  ,  index  on  D  }  cutoff  values  of  6  }  index  on  A 
Outputs  -  R  ,  T  ,  jZf  ,  t 


-  3  - 
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D.2  R  and  F  for  3  <  V  <  00 

Inputs  -  index  on  v  ,  index  on  D  ,  cutoff  values  of  9  ,  index  on  A 
Outputs  -  R  ,  T  ,  J2f  ,  t 

D.3  R  and  T  for  v  =  °° 

Inputs  -  index  on  D  ,  cutoff  values  of  B  9  index  on  A 
Outputs  -  R  ,  V  ,  0  ,  t 

D. 4  Subprogram  for  Printing  R  and  T 

Inputs  -  R  ,  T  ,  JZf  ,  ♦ 

Outputs  -  R  ,  T  ,  0  ,  f  tables  in  final  form 

E.  Programs  for  Computing  A  and  R  Cutoffs 

E .  1  A  and  R  for  v  <  00 

c  c 

Inputs  -  index  on  v  ,  index  on  D 
Outputs  -  6c  ,  Aq  ,  R^ 

E.2  A  and  R  for  v  =  00 

c  c 

Inputs  -  index  on  D 
Outputs  -  6c  ,  Aq  ,  R^ 

F.  Subprogram  for  Aitken1 s  Interpolation  (modification  of  SHARE  Distribution 
Number  355*^08) 

Inputs  -  table  of  arguments  and  corresponding  function,  new  argument 
Outputs  -  new  function 


G.  Supplementary  Programs  for  Computing  Optimal  Sampling  Functions 


G.1  Optimal  Sample  Size  -  First  Search  Procedure 

Inputs  -  v  ,  y  D  ,  A  ,  n  (maximum  number  of  trials) 

Outputs  -  R 


G.2  Optimal  Sample  Size  -  Second  Search  Procedure 

Inputs  -  v  ,  B(^,-|-v)  j  D  j  A  j  n  (maximum  number  of  trials) 

Outputs  -  R 

G.3  Optimal  Sample  Size  -  Newton-Raphson  Procedure 
Inputs  -  V  ,  B('lb'iv)  ,  D  ,  A  ,  trial  R 

Outputs  -  R 
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G.4  Net  Gain  for  Given  Optimal  Sample  Size,  Real  v  <  «»  }  t  <  v 

Inputs  -  v  ,  ,  D  ,  A  ,  R 

Outputs  -  r 

2 

G.5  Net  Gain  for  Given  Optimal  Sample  Size,  Real  v  <  00  ,  t  >  v 

Inputs  -  v  3  ,  D  ,  A  ,  R 

Outputs  -  r 

k.  Computer  Programs 

Computer  programs  summarized  in  Section  3  are  contained  in  this  section. 
Each  program  is  identified  by  title,  and  the  pages  are  numbered  within  programs. 
Comment  cards  are  liberally  used,  particularly  in  the  larger  7^90  programs  for 
computing  the  Student  distributions  and  the  main  Bayesian  tables. 
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Program  A. 1 ,  page  1 

A .  1  Standardized  Student  Functions  for  y  <  ea 

C - FS.  GS*  AND  FL5 »  NU  LESS  THAN  INF  »  T-0(. 01)10 

C 

PRINT  5 

5  FORMAT ( 3 1HOMOUNT  OUTPUT  TAPE  ON  LOGICAL  8) 

PAUSE 

C 

D  DIMENSION  FLNUT ( 19 ) »STVT ( 19 > .BETAT < 19 ) 

D  DIMENSION  T ( 1001 ) .FS( 1001 ) »GS( 1001 ) »FLS< 1001 ) 

D  DIMENSION  TSI91 > »FSS(91) 

C 

D  FLNUT ( i >*1* 

D  FLNUT (2 )“1*5 

D  FLNUT<3>*2. 

D  FLNUT 14 ) *3* 

D  FLNUT I  5 )*4. 

D  FLNUT { 6 ) *5. 

D  FLNUT ( 7 ) *6. 

D  FLNUT ( 8 ) *7. 

D  FLNUT { 9 )*8. 

D  FLNUTilOI«9. 

D  FLNUT ( 11 )*10« 

D  FLNUT ( 12  >*12  * 

D  FLNUT ( 13 )=15  • 

D  FLNUT ( 14)*20* 

D  FLNUT! 155*24. 

D  FLNUT 1 16) *30* 

D  FLNUT C 17) *40* 

D  FLNUT 118) *60* 

D  FLNUTi 191*120. 

C 

D  STVT(15*0. 31725  51743  05535  69877  E-Ol 

D  STVT 12) *0.11829  67755  68107  78653  E-Ol 

D  STVT135*  0. 492622 851 1662845 E-02 

D  STVTI4) *  0.10641995 292 07075 E-02 

D  STVT55)*  0.28100181 135 79956 E-03 

D  STVT 16 ) =  0.8547378787148 180E-04 

D  STVT <7  5*  0. 289 5991 3774768 13E-04 

D  STVT  ( 8 )  *  0. 10697 10 1445  3  8641 E-04 

D  STVT (9 )  =  0* 42440907638 14246 E-05 

D  STVT ( 10 ) *0. 17891187 15962368E-05 

D  STVT'  U)*0*7947765877982060E-06 

D  STVT (12) *0. 1790661843836 187E-06 

D  STVT ( 13 ) *0* 2 498449 071 462472 E-07 

D  STVT ( 14 ) *0.1581890879357 194E-08 

D  STVT  ( 15  !•  *0*  2  457762  7638  37653  E-09 

D  STVT ( 16  >*0.22876257041 14807E-10 

D  STVT (17 ) *0 • 9656558502 05 7788 E- 12 

D  STVT l 18 ) *0. 1068842957337041E-13 

D  STVT  ? 19 ) *0 » 8 56972 058768050 8 E- 17 

C 

D  SETATil)*  0.3141592653589793E+01 

D  BETA?  C 2  3  *  0.2396280469471184E+01 

D  8ETATI3)*  C.2000000000000000E+01 

D  BE TAT <4 )*  0. 1570796326794897E+01 

D  BETATI5)*  0. 1333333333333333E+01 

D  BETAT<6>*  0.1178097245096172E+01 

D  BETAT ( 7) *  0. 1066666666666667E+01 
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BETAT(8)*  0 • 98 1 7477042468 104E+00 
BETAT ( 9 ) *  0*9 142857142857 143 E+00 
BETAT  4X0) *0» 85902924 121 5959 1E+00 
BETAT (11 ) *0 ,8126984 12698412 7E+00 
BETAT ( 12 ) *0  * 7388 167388 167388E+00 
BETAT< 13 )*0.6580777580029401E+00 
BETAT  C 14 )»0» 567 5463855030419 E+00 
BETAT ( 15 )  *0.5 17019481 6 176779E+00 
BETAT ( 16 )*0*46l474553400974lE+00 
BETAT ( 17 )  *0. 3988 1 730689488 10E+00 
BETAT ( 18 )*0*3249554203948302E+00 
BETAT  < 19 ) *0 .2293000137934539E+00 

IT*  0 
T ( 1 ) *1 T 

DO  10  I  *  2*1001 
IT*IT  +1 
TT*IT 

10  Till  *  TT/100. 

WRITE  OUTPUT  TAPE  6*80* ( T< I) ,1*1*1001 ) 

IT*  9910 
TT*IT 

TSI 1 )*  TT/1000. 

DO  15  I*  2*91 
IT  *  IT  +  1 
TT-  IT 

15  TSCI)  *  TT/1000. 

WRITE  OUTPUT  TAPE  6*80*  (  TS(  I.)  »  1*1.91 ) 

-  READ  NUMBER  OF  NUS 

READ  INPUT  TAPE  5.20.NUMB 
20  FORMAT (1 10) 

-  MAIN  INDEX 

DO  500  IMAIN  *1 *NUMB 

-  READ  NU  INDEX  TO  SPECIFY  NU 

READ  INPUT  TAPE  5.25.INU 
25  FORMAT! 110) 

FLNU*FLNUT(INU) 

STV*STVT  <  INU) 

BETA*BETAT ( INU ) 

GS(1001)  *  STV 

-  COMPUTE  FS 

TERM1*  l./(BETA*SQRTF(FLNU) ) 

EXPr  ~( ( FLNU+1. ) /2. ) 
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DO  30  1*1*1001 

030  FS  Cl)*TERMl*  Ul.-MT  ( I  )**2  > /FLNU) **EXP ) 

C 

C  -  COMPUTE  GS 

C 

DO  35  1*1*91 

035  FSS< I)*TERM1*  (  ( l.  +  <TS(n**2)/FLNU>**EXP) 

C 

DO  50  1*2*10 

J«  11-1 
X*  10*IJ-1) 

0.  Tl*16067  •  *<FSS(K+mFSS(K+lln 

0  T 2*1063 00*  ‘  *( FSS<  K+2 1+FSSl K+10)  ) 

D  T3*48525*  »(FSS(K+3)4FSS(K+91> 

D  T 4*272400*  *<FSS(K+4)4FSS(K+8 >1 

D  T5*260550*  *(FSStK+5)+FSS<K+7) 1 

0  T 6*427366  *  *  C  FSS ( K+6 ) ) 

0  ATI*  C 5* *0*001 ) 7299376* 

D  AT2*Tl+T2-T3+T  4-T5+T6 

D  A*AT1*AT2 

L*1002-l 
M*1002-l4l 
050  G5(L)*GS(M)+A 

C 

DO  60  1*11*1001 
J*  1001-1 

0  Tl*16067*  #{FS  (J+ll+FS  (J+ll)) 

D  T2*106300.  *( FS  (J+2J+FS  (J+10M 

D  T3*48525.  *(FS  (J+31+FS  (J+9>) 

0  T4*272400.  *(F$  U+4J+FS  <J+8>) 

D  T5*260550.  *(FS  <J+5)+FS  <J+7)> 

D  T6*427368.  *(FS  (J+6)) 

D  ATI* ( 5**0*01  1/299376. 

0  AT2*Tl+T2-T3+T  4-T5+T6 

D  A-AT1*AT2 

K*1002-l 
L-1002-I+10 
060  GS(K )*GS(L )+A 

C 

C  -  COMPUTE  FLS 

C 

D  BOT  *  FLNU-1. 

DO  70  1*1*1001 

D  70  FLS1 1 )* ( (FLNU+  T ( I ) **2* ) /BOT >*FS<I)  -  T(I)*  GS«I> 
C 

C  -  WRITE  FS*GS*FL5 

C 

WRITE  OUTPUT  TAPE  8*  80*(FS< l ) .1*1.1001 > 

WRITE  OUTPUT  TAPE  8.  8o * (GSC 1 ) ♦ 1 *1 .1001 ) 

WRITE  OUTPUT  TAPE  8.  80 *(FLS( I> • 1*1 *1001 ) 

80  FORMAT  <5E20*8) 

C 

C - INDEX 

C 

500  CONTINUE 
C 

C - END 
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C 

END  FILE  8 
PRINT  999 

999  FORMAT ( 20H1  DISMOUNT  LOGICAL  8) 
PAUSE 
CALL  EXIT 
END 


-  q  - 
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A. 2  Standardized  Student  Functions  for,  v  =  » 


Program  A. 2,  page  > 


C  -  FS.  GS»  ANr  FLS*  NU  EQUALS  INF (NORMAL) .  T*0( .01)10 

C 

PRINT  5 

5  FORMAT ( 31HOMOUNT  OUTPUT  TAPE  ON  LOGICAL  8) 

PAUSE 

C 

D  DIMENSION  FLNUT ( 19 ) »STVT ( 19 ) .BETAT ( 19) 

D  DIMENSION  T ( 1001 ) »FS( 1001 ) »GS< 1001 ) »FLS< 1001 ) 

D  DIMENSION  TS ( 91 ) »FSS( 91 ) 

C 

IT*  0 

D  T(1)*IT 

DO  10  I  *  2*1001 
IT=IT  +  1 
D  TT= IT 

D  10  Til)  «  TT/100. 

C 

WRITE  OUTPUT  TAPE  6*80* I T ( I )• 1*1,1001 ) 

C 

IT*  9910 
D  TT*IT 

D  TS( 1 )*  TT/1000. 

DO  15  I*  2*91 
IT  *  IT  +  1 
D  TT*  IT 

D  15  TS< I )  =  TT/1000* 

C 

WRITE  OUTPUT  TAPE  6*80 » ( TS( I ) *  1*1 *91 ) 

C 

C  -  COMPUTE  FS 

C 

D  C0N=1./SQRTF<2. *3. 141592653589793) 

C 

DO  30  1*1.1001 

D30  FS!  I )*CON  *  EXPF (  (-.5)  *  T(I)**2) 

C 

C  -  COMPUTE  GS 

C 

GS(lOOl)*. 76198  53024  16059  16  E-23 
C 

DO  35  1*1.91 

D35  FSSt  I ) *CON  *  EXPF!  I-.5)  *  TSU>**2> 

C 

DO  50  1*2.10 

J*  11-1 
K*  10*( J— 1 ) 

D  T  1*16067 .  *(FSSOC+l)+FSS(K+llh 

D  T2*106300.  *(FSS(K+2)+FSS(K:+10)  ) 

D  T3*48525 •  *( FSS( K+3 )+FSS( K+9 ) ) 

D  T 4*272400.  * ( FSS! X+4 )+FSS( K+8 ) ) 

D  T5*260550»  *(FSS(K+5)+FSS( K+7) ) 

D  T  6*427368 •  *(FSS(K+6)) 

D  ATl*(5.#0.001)/299376. 

D  AT2*T1  +  T2-T3+T  4— T5+T6 

D  A*AT1*AT2 

L-1002— I 
M=1002-I+l 
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D50  GS(L)=GS(M)+A 

C 

DO  60  1*11*1001 
J=  1001-1 

D  Tl=16067.  * ( FS  (J+1)+FS  (J+ll)) 

D  T2*106300.  * ( FS  (J+2J+FS  (J+10U 

D  T3=48525«  *<FS  (J+3)+FS  (J+9H 

D  T4=272400.  *(FS  ( J+4J+FS  (J+8)> 

D  T5=260550*  *(FS  <J+5)+Fs  <J+7)J 

D  T6=427368*  *<FS  ( J+6 ) ) 

D  AT1=(5.*0.01  )/299376. 

D  .  AT2*T1+T2-T3+T4-T5+T6 

D  A=AT1*AT2 

K=1002-I 
L=1002-I +10 
D60  GS<K)=GS(L>+A 

-  COMPUTE  FLS 

DO  70  1=1*1001 

0  FLS(I)=FSU)  -  T  ( I  )  *GS(  I ) 

— -  WRITE  FS »GS*FLS 

WRITE  OUTPUT  TAPE  8*  80 » < FS< I ) . 1=1 » 1001 ) 
WRITE  OUTPUT  TAPE  8,  80 » (GS( I) » I *1*1001 ) 
WRITE  OUTPUT  TAPE  8*  80 * < FLS ( I) » I *1 *1001 ) 
80  FORMAT  (5E20.81 

- END 

END  FILE  8 
PRINT  999 

999  FORMAT ( 20H1  DISMOUNT  LOGICAL  8) 

PAUSE 
CALL  EXIT 
END 


11 
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B.  Program  for  Printing  Tables 


^s*  ’  ,Gs*  ancl  Ls* 


-  F5.  GS»  AND  FL5  EDIT  AND  PRINT 

PRINT  10 

10  FORMAT ( 30H1M0UNT  INPUT  TAPE  ON  LOGICAL  3) 

PRINT  15 

15  FORMAT ( 31H1MOUNT  OUTPUT  TAPE  ON  LOGICAL  8) 

PAUSE 

DIMENSION  DISTU00I) 

DIMENSION  T( 11 ) »TT  1 10 ) 

DIMENSION  Pl(ll)  ».P2  ( 1 1 )  »P3 ! 11 )  »DIF ( 5 ) 

DIMENSION  IP3 ( 11 ) * ID( 11 ) 

DIMENSION  I P4( 1 1 > 

105  FORMAT (3110) 

111  FORMAT ( 1H112X» 1HF 1 3 *4 IX » 16HSTUDENT  0RDINATE41X * 1HFI 3 ) 

112  FORMAT ( 1H112X * 1HGI 3 *40X» 18HSTUDENT  RIGHT  TAIL40X.4X) 

113  FORMAT! 1H112X.4X*  4QX* 19HSTUDENT  LINEAR  LOSS39X *  1HL 13 ) 

125  FORMAT!  18X*  78HEXP  0  1234  5 

1  6  7  8  9  10) 

130  FORMAT! 5E20. 8) 

131  FORMAT ( £20*8  »4X ) 

139  FORMAT! 10X*6H  O.OI5*4X»I4»4X*I4*2X*I4*2X*I4*2X*I4»4X*I4*4X* I4*2X 

1*I4»2X»I4»2X*I4»4X»I4»2X*5I4) 

140  FORMAT ! 10X*F6« 1 ♦  15 »4X* 14 *4X f I4*2X. I4»2X* I4*2X* I4»4X» I4t4X» I4*2X 

1* I4»2X» I4«2X* I4»4X* I4»2X»5I4) 

141  FORMAT! 10X*F6.1*  15 »4X » 14 *4X» l4t2X » I4.2X. I4»2X. I4.4X. I4.34X* . 

1  2X  *5 1 4 ) 

142  FORMAT! 10X.6X*  I5*4X*34X*  I4»4X*I4»2X 

1*I4»2X*I4*2X*I4»4X*I4*2X*5I4) 

145  FORMAT! 1H  ) 

READ  INPUT  TAPE  5*105*NUMNUS 
WRITE  OUTPUT  TAPE  6.105 ,NUMNUS 

DO  2000  IMAINA=1 .NUMNUS 

-  ROUTINE  FOR  EACH  NU 

READ  INPUT  TAPE  5 *  105 .NU .LINSW »L INTOT 
WRITE  OUTPUT  TAPE  6 . 105 .NU. LINSW .L I NTOT 
C 

DO  1800  IMAIN8=1 .3 


C 

IFCIMAINB  -2)211*212*213 


211 

WRITE 

OUTPUT 

TAPE 

8*111 *NU*NU 

WRITE 

OUTPUT 

TAPE 

8*145 

WRITE 

OUTPUT 

TAPE 

8.125 

WRITE 

OUTPUT 

TAPE 

8*145 

GO  TO 

230 

212 

WRITE 

OUTPUT 

TAPE 

8.112. NU.NU 

WRITE 

OUTPUT 

TAPE 

8.145 

WRITE 

OUTPUT 

TAPE 

8*125 

WRITE 

OUTPUT 

TAPE 

8.145 

GO  TO 

230 

213 

WRITE 

OUTPUT 

TAPE 

8. 113. NU.NU 

WRITE 

OUTPUT 

TAPE 

8,145 

12 
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WRITE  OUTPUT  TAPE  8.125 
WRITE  OUTPUT  TAPE  8.145 
GO  TO  230 

-  ROUTINE  FOR  EACH  DISTRIBUTION 

230  READ  INPUT  TAPE  3.130. CDISTC I ) .1*1.1001) 
WRITE  OUTPUT  TAPE  6 .105 ♦ IMAINA .NU. IHAINB 

ARG=0.0 

NL*1 

BN=10.**30 
I LINE*1 

DO  1000  IMAINC*1.80 

DO  240  JPR=1 .11 
KPR*10*(  IMAINC— 1 )+JPR 
PHJPR)*DIST  (KPR) 

240  CONTINUE 

-  DECISION  1  OR  2  LINES 

1F( ILINE-LINSW) 300.600.600 

-  1  LINE 

300  ILINE*ILINE+1 
A1=P1 ( 1 ) *BN 

A2=  .43429448  *  LOGF(Al) 

1 1=A2" 

IE=MIl-30)  +1 
I E2=  XABSF  (IE) 

IF( IE-0  1330.320.310 

310  FMULT  =  10.**(-IE2) 

DO  315  LPRa 1 .11 
315  P2ILPR) «  PI ( LPR ) *  FMULT 
GO  TO  340 

320  DO  325  LPR=>1.11 
325  P2(LPR) =  Pl(LPR) 

IE*+0 
GO  TO  340 

330  FMULT  *10. **  IE2 
DO  335  LPR-1.11 
335  P2 (LPR) =P1 ( LPR ) *FMULT 
GO  TO  340 

340  DO  345  J=l*ll 

345  P3  ( J ) *P2 ( J)  *  10000. 

DO  360  J«l.ll 
AX=P3 ( J ) 

1  AY*P3 ( J ) +1. 

I  A*AY 
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AZ=  I A 

IF< (AZ-AX>-(AY-AZ> >350.350*355 
350  I P3 ( J ) = AY 
GO  TO  360 
355  I P3 ( J ) = AX 
GO  TO  360 
360  CONTINUE 

DIFT  =0. 

DO  365  J»2»ll 

365  DIFT=DIFT  +(  P3 U > -P3 < J-l ) ) 

DIFM*ABSF(DIFT/10.) 

D I F ( 1 )=DIFM*»1 
DIF ( 2 ) sDIFM*«2 
DIF(3)*DIFM*.3 
DIFJ4)=DIFM*.4 
DIF  (  5 )  *>DIFM**5 
DO  380  J«1 .5 
AX*DIF ( J ) 

AY=DIF(J>+1. 

I  A=AY 
AZ=IA 

I F ( ( AZ-AX ) - ( AY-AZ 1 >370.370.375 
370  ID( J )=AY 
GO  TO  380 
375  ID( J )*AX 
GO  TO  380 
380  CONTINUE 

DO  999  J*l»ll 
999  IP4(J)=IP3{J) +10000 

IF  (ARG-0 *0)381. 381  *382 

381  WRITE  OUTPUT  TAPE  8.139.IE.  ( IP4< J) » J=1 . 11 ) , ( ID< J > . J*1 .5 > 

GO  TO  384 

382  WRITE  OUTPUT  TAPE  8.  140.  ARG.  I E.  C IP4I  J )  .  J*1 . 11 )  .  (IOU )  .  J«1 .5  ) 
GO  TO  384 

384  NL=NL+1 

IF (NL-6  >385*390.390 

385  GO  TO  900 

390  WRITE  OUTPUT  TAPE  8.  145 
NL=  1 

GO  TO  900 

-  2  LINE 

600  ILINE=I LINE+2 
1 2L  =  1 
M=1 
MP*2 
N=6 

GO  TO  630 
C 

610  I 2L*2 
M=6 
MP*7 
N  =  1 1 
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GO  TO  630 
C 

630  A1=P1(M)*BN 

A2=.43429448*LOGF ( A1 ) 

1 1  =  A2 

I £=  (  I 1-30 )  +  l 
I  E2*XABSF (IE) 

C 

I F ( I E—0  >680*670*660 
C 

660  FMULT=10.**(-IE2) 

DO  661  J=M»N 

661  P2(J)*P1(J)*FMULT 
GO  TO  700 

C 

670  DO  671  J=M.N 

671  P2(J)=P1(J) 

IE=+0 

GC  TO  700 
C 

680  FMULT*10.**IE2 
DO  681  J=M»N 

681  P2(J)=P1(J)*FMULT 
GO  TO  700 

C 

700  DO  710  J«M.N 
710  P3 { J ) *P2 1 J ) *10000* 

DO  740  J=M*N 
AX=P3( J) 

AY=P3< J)+l. 

I  A=AY 
A2=  I A 

1 F ( { AZ-AX )-<  AY-AZ ) >720*720*730 
720  I P3 ( J ) =AY 
GO  TO  740 
730  I P3 ( J ) =AX 
GO ‘TO  740 
740  CONTINUE 
C 

D IFT=0* 

DO  750  J=MP*N 

750  DIFT=DIFT+(P3< J)-P3( J-l ) > 

D IFM=ABSF ( DIPT /5» ) 
DIF(l)*DIFM*.l 
DIF(2)=DIFM*.2 
D IF ( 3 ) =DIFM*«3 
DIF(4)=DIFM«.4 
DIF ( 5 ) =DI FM*«  5 
DO  790  J=1 *  5 
AX=DIF( J) 

AY=DIF( J)+l. 

I  A=AY 
AZ=  I A 

IFC  lAZ-AX)-(AY-AZ) >760*760*770 
760  ID(J)*AY 
GO  TO  790 
770  1D( J )=AX 
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GO  TO  790 
790  CONTINUE 
C 

GO  TO  (801.806) .I2L 
C 

801  DO  802  J=1  *6 

802  IP4(J)=IP3(J) +10000 

WRITE  OUTPUT  TAPE  8 . 141 .ARG , IE. ( I P4 ( J) • J*1 .6 ) » ( ID ( J ) . J=1 .5 ) 
GO  TO  610 
C 

806  DO  807  J=6 .11 

807  IP4( J)=IP3( J)+10000 

WRITE  OUTPUT  TAPE  8.142.  IE.  ( IP4 ( J) . J=6. 1 1 ) . ( ID( J ) . J= 1  .5 ) 

NL=NL+1 

IF(NL- 6)810. 820. 820 
810  GO  TO  900 

820  WRITE  OUTPUT  TAPE  8.145 
NL  =  1 

GO  TO  900 

-  DECISION  ON  END  OF  TABLE-  INDEX  ON  LINES 

900  IF(ILINE-LINTOT)920. 950*950 
920  ARG=ARG+0.1 
GO  TO  1000 
950  ARG=ARG+0.1 
GO  TO  1050 

1000  CONTINUE 

-  LAST-VALUE  IN  DISTRIBUTION 

1050  A1=P1 ( 1 1 )*BN 

A2=  *43429448  *  LOGF(Al) 

11=  A2 

I E- 1 1 1—30)  +  1 
IE2=  XABSF(IE) 

C 

IF ( IE-O) 1070.1060.1070 
1060  I£=0 

PL1=P1 ( 11 ) 

GO  TO  1080 

1070  PL1=P1 (11)  *  (10.**  I E2 ) *  10000. 

GO  TO  1080 
C 

1080  AX=PL1 

AY=PL1  +1. 

IA=AY 
AZ  =  I  A 

IF  ( (AZ-AX)-(AY-AZ) )  1085.1085.1090 
1085  I PL=AY 

GO  TO  1100 
1090  IPL=AX 

GO  TO  1100 
C 

1100  WRITE  OUTPUT  TAPE  8.  140*  ARG.IE.IPL 
C 
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-  MAIN  INDEXES 

1800  CONTINUE 
2000  CONTINUE 
C 

END  FILE  8 
PRINT  3000 

3000  FORMAT! 2 5H1 DISMOUNT  LOGICAL  3  AND  8) 
PAUSE 
CALL  EXIT 
END 
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-  SERIES  EXPANSION  FOR  B!1/2,1/2NU> 

DIMENSION  RT!2)»IR!2)»GRT<2) ,GR<2) 
C 

5  FORMAT (110) 

10  FORMAT (E40.20,I10) 

20  FORMAT! E40. 20 ) 

25  FORMAT! 1H  ) 

C 

SQPI=SQRTF( 3.1415926535897932385 ) 
PRINT  20.SQPI 
C 

READ  5»NUMB  ' 

DO  1000  1MAIN=1 *NUMB 
C 

READ  10,FLNU,NZ 
PRINT  10.FLNU.NZ 
X=«5*FLNU 
C 

IX=X 
TIX=IX 
r=X-T I X 

IF! R- .5)1 30 *140.140 
130  RT!1)=R+1. 

I R ( 1 )=IX-1 
RT 1 2 )=R+.5 
I R ( 2 ) =1 X 
GO  TO  200 
140  RT ( 1 ) =R 
IR(1)*IX 
RT  ( 2')=R+.5 
I R  (  2  )  =  I X 
GO  TO  200 
C 

200  DO  400  1=1,2 

I F ( RT ! I )-l. >250,225,250 
225  GRT ( I ) = 1 • 

GO  TO  400 
250  ZT =0  ■ 

Zl  =  l. 

DO  300  J=1 » NZ 
F  J=J 

READ  INPUT  TAPE  2,10*ZETA 
Zl=( l.-RT ! I ) )»Z1 
Z=(Z1*ZETA)/FJ 
300  ZT=ZT+Z 

GRT! I ) =EXPF ! ZT ) 

REWIND  2 
400  CONTINUE 

PRINT  20, GRT! 1) »GRT<2) 

C 

DO  600  1=1,2 
GR  < I )*GRT (  I ) 

RMULT=RT( I ) 

C 

1 RT= IR 1 1 ) 

DO  500  J=1 *  I RT 
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GR(I)=GR(I)*RMULT 
RMULT«=RMULT  +  1. 

500  CONTINUE 
C 

600  CONTINUE 

BETA=(SQPI*GR( 1 ) )/GR(2) 

PRINT  20 »GR ( 1 ) »GR( 2  I *BET A 
PRINT  25 
C 

1000  CONTINUE 
END 
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C.2  Finite  Series  Expansion  for  ,  Integral  Even  v  ■  *• 


-  SERIES  EXPANSION  FOR  GS»  INTEGRAL  EVEN  NU 

-  FS»  GS»  AND  FLS 

10  FORMAT  < 1 10 ) 

17  FORMAT ( I10»E30«20»F10«2) 

18  FORMAT ( E30.20 ) 

19  FORMAT ( 1H  ) 

READ  10»NUMNUS 
PRINT10  »NUMNUS 
DO  450  I NU=  1 *NUMNUS 

READ  17»NUfBETA,T 
PR INTI 7  »NU*BETA*T 
FLNU=NU 

-  DENSITY  FS 

Z=(1./SQRTF(FLNU) >*<1./BETA) 
EXP=(-.5)*(FLNU+1. ) 

FS=Z*< l.+T**2/FLNU)**EXP 
PRINT  18  »FS 

-  RIGHT  TAIL  GS 

IMAX=NU/2 

CONA=FLNU/ ( FLNU+T**2 ) 

TERMU=1 • 

SUMU=1. 

DO' 440  1=2  *  I MAX 
FI  =  I 

TINCR=(2.*FI-3.)/(2.*FI-2.> 

TERMU=T I NCR*CONA*TERMU 
SUMU=SUMU+TERMU 
440  CONTINUE 

CONB=T /SQRTF ( FLNU+T**2 ) 

GS=.5*( l.-CONB*SUMU) 

PRINT  18  *GS 

-  LINEAR  LOSS  FLS 

FLS= ( ( FLNU+T **2  > / ( FLNU~1 • ) ) *FS  -  T*GS 
PRINT  18.FLS 
PRINT  19 

450  CONTINUE 
END 
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C.3  Finite  Series  Expansion  for  Gg^(tjv)  ,  Integral  Odd  v  <  °° 


-  SERIES  EXPANSION  FOR  GS.  INTEGRAL  000  NU 

-  FS»  GS,  AND  FLS 

17  FORMAT! I 10,E30.20»F10.2 ) 

18  FORMAT (E30. 20) 

21  READ  17  »NU»BETA*T 
PRINT17»NU,BETA*T 
FLNU*NU 

-  DENSITY  FS 

Z* ( 1 . /SQRTF ( FLNU ) ) » ( 1 */BETA ) 
EXP*C-.5)*(FLNU+1.) 
F5*Z*(1.+T**2/FLNU)**EXP 
PRINT  18 ,FS 

-  RIGHT  TAIL  GS 

IMAX*(NU-l)/2 
C0NA=FLNU/CFLNU+T**2)  ' 

IF (NU— 3 ) 41 0,41 5 *420 
410  SUMV “0  » 

GO  TO  490 
415  SUMV*C0NA 
GO  TO  490 
420  TERM*C0NA 
SUMV=C0NA 

DO  480  I *2  *  I MAX 
FI*  I 

TINCR=(2.*FI-2.)/(2.*FI-l.) 

T  ERM=T I NCR*CONA*TERM 
SUMV=SUMV+TERM 
480  CONTINUE 

-  ARCTAN  OF  SQRTF ( FLNU/T ) 

490  X*SQRTF ( FLNU ) /T 

PI *3* 1415926535897932385 
C 

I F(X— 1*5)491 ,491 ,496 
C 

491  C0Nl=X**2/( l.+X**2) 

TERM=1. 

SUM*1 * 

SUMB4*1. 

C 

DO  493  1*2,100 
F I  *  I 

TINCR* ( 2»*FI-2  *)/(2**FI-l*) 

TERM=T 1 NCR*C0N1#TERM 

SUM*SUM+TERM 

I F ( SUM-SUMB4 ) 494 , 494 ,492 

492  SUMB4=SUM 

493  CONTINUE 

494  C0N2«X/(1.+X*»2) 
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AT  ANX=CON2*SUM 
GO  TO  510 
C 

496  C0N1=1» /X**2 
TERM=1./X 
SUM=1*/X 
SUMB4=0. 

C 

DO  498  1=2  .100 
F I  =  I 

TINCR=-<<2.*Fl-3.)/<2.*FI-l.) > 
TERM=TINCR*C0N1*TERM 
SUM=SUM+TE'RM 
I F ( SUM-SUMB4 )  497 , 499 .497 

497  SUMB4=SUM 

498  CONTINUE 

499  ATANX=PI /2.-SUM 
GO  TO  510 

C 

510  PRINT  18.X.ATANX 

CONB=T/SQRTF (FLNU ) 

GS= < 1 . /P I)* ( AT ANX-CONB*SUMV ) 

PRINT  18 .GS 

-  LINEAR  LOSS  FLS 

FLS=( (FLNU+T**2)/(FLNU-1.) >»FS  -  T.*GS 
PRINT  18. FLS 
C 

GO  TO  21 
END 


22 
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,  t2  <  v 


-  SERIES  EXPANSION  FOR  GS»  NU  REAL*  CONVERGES  FOR  T**2  LESS  THAN  NU 

-  FS»  GS*  AND  FLS 

10  FORMAT (1 10) 

17  FORMAT ( 2E30*20»F10*4) 

18  FORMAT!  E30.20»I5) 

19  FORMAT ( 1H  ) 

READ  10,  NUMNUS 
PRINT10,  NUMNUS 
DO  500  I NU=1 » NUMNUS 

READ  17,TNU,TBETA*T 
PRINT17»TNU»TBETA,T 

-  DENSITY  FS 

2=(1./SQRTF(TNU) )#(1./TBETA) 

EXP«(-.5)*(TNU+1.) 

FS=Z*( l.+T**2/TNU)**EXP 
PRINT  17  *FS 

-  RIGHT  TAIL  GS 

CONA=T**2/TNU 
TERM*1 • 

SUM=1. 

SUMB4-0 • 

DO  340  J=2*500 
F  J=J  " 

T 1 NCR*- ( ( 2 . *F J-3 • ) * ( TNU+2. *F J-3 . ) ) / ( ( 2 • *F J-l . ) * ( 2 • *F J-2 . )) 

TERM*T I NCR*CONA*TERM 
SUM=SUM+TERM 
I F ( SUM-SUMB4 ) 338*345,338 
338  SUMB4=SUM 
340  CONTINUE 
GO  TO  999 

345  CONB= ( T /SQRTF ( TNU ) ) /TBET  A 
GS=.5-CONB#SUM 

PRINT  18  »GS* J 

—  LINEAR  LOSS  FLS 

FLS= ( ( TNU+T**2 ) / ( TNU—1 • ) )*FS  -  T*GS 
PRINT  17  *FLS 
PRINT  19 
C 

999  CONTINUE 
500  CONTINUE 
END 
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i  2 

C.5  Infinite  Series  Expansion  for  Gg^(t|v)  ,  Real  V  <  »  }  t  '•  v 


-  SERIES  EXPANSION  FOR  GS.  NU  REAL*  CONVERGES  AS  POWER  SERIES  IN 

- NU/ ( NU+T **2 )  FS*  GS*  AND  FLS 

10  FORMAT (1 10) 

17  FORMAT (2E30.20»F10.3) 

18  FORMAT (  £30.20,15) 

19  FORMAT ( 1H  > 

READ  10,  NUMNUS 
PRINT10 »  NUMNUS 
DO  500  INU*1 *NUMNUS 

READ  17*TNU*TBETA.T 
PRINT17*TNU*TBETA.T 

-  DENSITY  FS 

Z*  ( 1 */SORTF( TNU ) )*( l./TBETA) 

EXP*(-.5)*(TNU+1.) 

FS=Z*(1.+T**2/TNU)**EXP 
PRINT  17  *FS 

-  RIGHT  TAIL  GS 

CONA=TNU/ ( TNU+T**2 ) 

TERM-1. /TNU 
SUM*  1 • /TNU 
SUMB4*0* 

DO  340  J=2  *500 
F  J=J 

TINCR»( (2.*FJ-3.)*<TNU+2.*FJ-4.) )/( < 2.*FJ-2. >*< TNU+2.*FJ-2. ) ) 
TERM*TINCR*CONA*TERM 
SUM=SUM+TERM 
I F(SUM-SUMB4) 345*345 *338 
338  SUMB4=SUM 
340  CONTINUE 
GO  TO  999 

345  CONB*( l./TBETA) *CONA**( #5*TNU> 

GS*CONB*SUM 
PRINT  18*GS*J 

-  LINEAR  LOSS  FLS 

FLS* ( ( TNU+T**2 ) / ( TNU-1* ) )*FS  -  T*GS 
PRINT  17, FLS 
PRINT  19 
C 

500  CONTINUE 
999  CONTINUE 
END 
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—  NORMAL  RIGHT  TAIL  USING  LAPLACE  SERIES.  T  GREATER  THAN  2 

1  FORMAT ( I I0.E40.20 ) 

2  FORMAT ( E40.20 ) 

C 

READ  1 .NUMB 
DO  200  IMAIN=1 .NUMB 
C 

READ  l.N.T 
AN*N 

Z*AN/(T+AN+1 • J 
NM1*N-1 

DO  100  I *1 »NM1 
AI*I 

100  Z=( AN-AI >/( T+Z ) 

PI »3 .1415926535897392 
F»(EXPF(-.5*T**2) )/SQRTF(2**PD 
G»F/(T+Z) 

PRINT  2.G 
GNBS«1.-2.*G 
PRINT  2.GNBS 
C 

200  CONTINUE 
END 
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C.7  Series  Expansion  for  Gg#(t|v  =  °°)  ,  t  <  2 


-  NORMAL  RIGHT  TAIL*  T  EQUAL  TO  OR  LESS  THAN  2 

10  FORMAT (I10*F10»2) 

20  FORMAT (E40. 20) 

C 

READ  10»NUMB 
PRINT10.NUMB 
DO  300  IMAIN*1*NUMB 
C 

READ  10* ITER*T 
PRINT10* ITER»T 
PI *3 •1415926535897932 
CON*l./(SQRTF(2**PI )) 

SUMB4*0« 

TERM-T 

SUM*T 

DO  100  I  *1 » I  TER 
AI*I 

TINCR»-<2**Al-l.)/(2#*AI»(2.*AI+l.n 
TERM»TINCR*TERM*T**2 
SUMB4*SUM 
SUM-SUM+TERM 
I F ( SUM-SUMB4 ) 1 00  *  1 5 0 . 1 00 
100  CONTINUE 
150  G  *  .5-(SUM*CON> 

PRINT  10»I 
PRINT  20»G 
C 

300  CONTINUE 
END 
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D.  1  R  and  I*  for  y  =  1-^  ,  2 

-  OPTIMAL  SAMPLE  SIZE  AND  NET  GAIN.  NU-3/2.2  • 

PRINT  5 
PAUSE 


DIMENSION  CL  (121)  .CLGRHO <41 » 121 )  .CLGNGUl  .121 ) 
COMMON  CL. CLGRHO »C LONG 
DIMENSION  PHIA5(121).PSIAS(121) 

COMMON  PHI AS .PS I AS 


DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 


FLNUT ( 19 ) *STVT ( 19 ) .BETAT ( 19) 

T(100l).  FS < 1001 )*GS (1001 )»FLS( 1001) 

TS(9l)*FSS(91 ) 

L  (  11 ) 

D(  41 ) .THCTT (  41 ) »THCT (  4l).FLAMCT(  41).FLCT<  41) 
ATH( 360) tBL( 121) »BLAM( 121 ) 

ALAM ( 360 ) .CLAM ( 121) »CTH( 121) 

CRHOT (121) »CRHQ( 121 ). I  TER ( 121 ) .CGPR ( 121 ) 

CRH0AS( 121 ) .RHORAT (121) 

CNG( 121 ) 


5  FORMAT (38H1MOUNT  OUTPUT  TAPES  ON  LOGICAL  2  AND  8) 

6  FORMAT (25H1D I SMOUNT  LOGICAL  2  AND  8) 

10  FORMAT ( 3110) 

15  FORMAT! 5E20. 8) 

16  FORMAT (5F20. 8) 

21  FORMAT (1H1) 

22  FORMAT! 1H  ) 

51  FORMAT! 2 1H0STUDENT  DENSITY.  NU«I4) 

52  FORMAT ( 24H0STUDENT  RIGHT  TAIL.  NU-I4) 

53  FORMAT! 25H0STUDENT  LINEAR  LOSS*  NU-I4) 


FORMAT! 110H0 

CL 

CLAM 

CTl 

CRHO 

CLGRHO 

CNG 

CLGNG) 

FORMAT! 70H0 

RHORAT ) 

CL 

CRHO 

CRHOAS 

71  FORMAT ( F10.2 .3E20.8 .F10.4.E20.8  *F10*4) 

72  FORMAT (F10.2.3E20.8.I10) 

73  FORMAT (F10.2.3E20.8) 

81  FORMAT! 3HOD-F6. 2) 


SPECIFY  FLNUT*  STVT.  AND  BETAT 


FLNUT ( 2 )-l • 5 
FLNUT(3)«2. 

STVT (2) =0.11829  67755  68107  78653  E-Ol 
STVT ! 3 ) =  0. 4926228 51 1662845 E-02 

BETAT ( 2 ) *  0. 239628 04694 7 1184E+01 
BETAT ( 3 ) *  0 • 2QOQOQOOOOOOQOOQE+Q 1 

-  GENERATE  T  AND  TS 

IT-0 
D  T ( 1 ) “IT 

DO  110  1-2.1001 
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IT*  IT  +  1 
TT=IT 

110  T(I)  *  TT/100. 

IT  =  9910 
TS(l)*  IT/1000 
DO  115  1=2*91 

IT-IT+l 
TT=  I T 

115  TS( I )=TT/1000. 

—  D  LIST  01.1)4 

ID  =  0 
D  ( 1 ) =  ID 
DO  130  1=2*41 
ID* iD+1 
TD=  ID 

130  D(I)=  TD/10. 

-  ATH  LIST  360  VALUES  *9999  9999  TO  *01 


ATHI 

1 ) =.9999 

9999 

ATH( 

2) *.9999 

9998 

ATHI 

3) *.9999 

9997 

ATHI 

4)  =.9999 

9996 

ATHI 

5) =.9999 

9995 

ATHI 

6 ) *.9999 

9994 

ATHI 

7 ) =.9999 

9993 

ATHI 

8)*. 9999 

9992 

ATHI 

9 >=.9999 

9991 

ATHI 

10) *.9999 

999 

ATHI 

ID*. 9999 

998 

ATHI 

12 >*.9999 

997 

ATHI 

13)  *.9999 

996 

ATHI 

14)*. 9999 

995 

ATHI 

15)=. 9999 

994 

ATHI 

16)*. 9999 

993 

ATHI 

17) **9999 

992 

ATHI 

18)*. 9999 

991 

ATHI 

19) *.9999 

99 

ATHI 

20) *.9999 

98 

ATHI 

21) *.9999 

97 

ATHI 

22)*. 9999 

96 

ATHI 

23)=. 9999 

95 

ATHI 

24) =.9999 

94 

ATHI 

25)=. 9999 

93 

ATHI 

26)=. 9999 

92 

ATHI 

27 ) *.9999 

91 

ATHI 

28)*. 9999 

9 

ATHI 

29)-. 9999 

8 

ATHI 

30) *.9999 

7 

ATHI 

31)=. 9999 

6 

ATHI 

32)*. 9999 

5 

ATHI 

33)*. 9999 

4 

ATHI 

34)*. 9999 

3 

ATHI 

35)*. 9999 

2 
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ATHi  36)=. 9999  1 
ATH(  37)=. 9999 
DO  133  J=38*236 
133  ATH(  J)=ATH(  J-D-.0001 
ATH( 237 ) *.979 
DO  134  J=238 *266 

134  ATH( J)»ATH( J-D-.001 
ATHC  267 ) =.94 

DO  135  J=268  *360 

135  ATH( J)=ATH( J-1J-.01 
C 

WRITE  OUTPUT  TAPE  6»15t<  ATH( J) »J* 1*360) 
WRITE  OUTPUT  TAPE  6*16* < ATH (J ) »J=1 *360) 

-  BL  AND  BLAM  L«-2.(.l)10 

LT*-20 
BL(l)*-2. 

BLAM ( 1 ) *.01 
DO  160  J=2  *  121 
LT=LT+1 

IF(LT  -  0)155*150*155 
150  BL(J)»  +  0 
GO  TO  160 
155  FLT  *  LT 

BL( J  )*  FLT/10. 

160  BLAM(J)  =  10.**  BL ( J) 

-  READ  NUMNU  *  MAIN  INDEX  ON  NU  -  IOFNU 

READ  INPUT  TAPE  5*10*NUMNU 
WRITE  OUTPUT  TAPE  6*10*NUMNU 

DO  3000  IOFNU  *  1*NUMNU 

-  READ  INU*  NU*  AND  MAXIT 

READ  INPUT  TAPE  5 » 10* INU*NU*MAXIT 
WRITE  OUTPUT  TAPE  6*10* INU»NU*MAXIT 

FLNU=FLNUT( INU) 

5TV=STVT  < INU ) 

BETA=BETAT( INU) 

WRITE  OUTPUT  TAPE  6*15  »FLNU*STV *BETA 

-  INITIALIZE  CLGRHO  AND  CLGNG 

DO  200  1*1*41 
DO  190  J*1 .121 

CLGRHOI I  * J )  *  0. 

CLGNGU.J)  *0. 

190  CONTINUE 
200  CONTINUE 

-  COMPUTE  FS 

TERM1*  l./(BETA*SQRTF (FLNU) ) 
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D  EXP=-( ( FLNU+1* )  /  2  •  ) 

C 

DO  230  1=1.1001 

D  230  FS  (  I ) =TERM1*  1(1.+<T  ( I )**2 ) /FLNU) **EXP  J 
C 

C  -  COMPUTE  GS 

C 

D  GS (1001)  =  STV 

C 

DO  235  1*1.91 

D  235  FSS( I )=TERM1*  ( ( 1 *+(TS ( I )**2 ) /FLNU) **EXP ) 

C 

DO  250  1*2*10 
J=  11-1 
K*  10*(J-1) 

D  Tl-16067  •  *(FSS(K+1)+FSS(K+11)) 

D  T2*106300.  *(FSS(K+2)+FSS(K+10>) 

D  T3=48525.  * ( FSSI K+3 )+FSS( K+9 ) ) 

D  T 4*272400*  * t FSS< K+4)+FSS( K+8 ) ) 

D  T5*260550.  * ( FSS1K+5 J+FSSl K+7 ) ) 

D  T6*427368*  *<FSS(K+6)> 

D  ATl*<5.*0.001)/299376. 

D  AT2*T1+T2-T3+T4-T5+T6 

D  A=AT1*AT2 

L*1002-I 
M=1002-I+l 
D  250  GSCL )*GS(M)+A 
C 

DO  260  1*11*1001 
J*  1001-1 

D  Tl  =  16067.  *<FS  (J+D+FS  (J+1U) 

D  T2=106300*  * (FS  IJ+2)+FS  (J+lO)) 

D  T3*48525.  *(FS  (J+3J+FS  ( J+9 )) 

D  T4=272400.  *(FS  (J+4J+FS  ( J+8 ) ) 

D  T5*26055Q*  *(FS  (J+5H-FS  (J+7)) 

D  T6*427368*  *(FS  (J+6>) 

D  ATI* ( 5**0*01  J/299376. 

D  AT2*T1+T2-T3+T4-T5+T6 

D  A*AT1*AT2 

X*1002-I 
L*1002-I+10 
D  260  GS(K)=GS(L)+A 
C 

C  -  COMPUTE  FLS 

C 

D  BOT  *  FLNU— 1 • 

DO  270  1=1*1001 

D  270  FL5( I J  * ( (FLNU+  T(I)**2. > /BOT)*FS( I )  -  T ( I )*  GS( I) 
C 

C  -  PRINT  STUDENT  DISTRIBUTIONS 

C 

WRITE  OUTPUT  TAPE  6»51*NU 

WRITE  OUTPUT  TAPE  6.15*(FS  ( I ) *1=1*1001*100) 

WRITE  OUTPUT  TAPE  6.52.NU 

WRITE  OUTPUT  TAPE  6.15. (GS  ( I) *  1*1 . 1001 .100 ) 

WRITE  OUTPUT  TAPE  6.53.NU 

WRITE  OUTPUT  TAPE  6 . 15 . ( FLS ( I ) ♦ 1*1.1001 .100 ) 
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-  READ  INDEX  IOFD  INDEX  ON  D 

READ  INPUT  TAPE  5. 10* I0FD1 , I0FD2. I0FD3 
WRITE  OUTPUT  TAPE  6 .10* I0FD1 ♦ IOFD2 » IOFD3 

DO  2000  IOFD  =  !OFDl.IOFD2*IOFD3 

-  READ  AND  PRINT  THCT(IOFD)*  IOFLl ♦ IOFL2 ♦ IOFL3 

READ  INPUT  TAPE  5*15»THCT< IOFD) 

READ  INPUT  TAPE  5* 10* IOFLl » IOFL2 . IOFL3 

WRITE  OUTPUT  TAPE  6*21 

WRITE  OUTPUT  TAPE  6*15»THCT( IOFD) 

WRITE  OUTPUT  TAPE  6*10* IOFLl ♦ IOFL2 ♦ IOFL3 

-  COMPUTE  ALAM  *  F(ATH)  WHERE  ATH  *  OR  GREATER  THAN  THCT(IOFD) 

-  NOTE  ATH  IS  STORED  IN  DECREASING  ORDER 

DO  400  J=1  *  360 

IFCATHI J>-  THCTI IOFD) ) 350 *350 *370 

350  JLONG= J-l 
GO  TO  425 

370  DTH  *  D ( IOFD ) /  ATH(J) 

L(l)=  XLOCFIL(l)) 

L (2 ) *  XLOCF ( T ( 1 ) ) 

L ( 3 ) *  XLOCF (FS( 1) ) 

L  ( 4 )  *  1 
L  ( 5 )  *  1 
L  ( 6 )  *  5 
L ( 7 ) *  1001 
L  t  8 )  *  0 

FSDTH  *  TABF ( DTH »L ( 1 ) ) 

I=L( 8) 

GO  TO( 390*375 ) *  I 

375  PRINT  380*  IOFD.  J*  DTH 

380  FORMAT ( 15H1TFS( D/ATH)  BAD2 1  5 »E20 .8 ) 

GO  TO  9999 

390  TERM1*(FLNU-1. )/<FLNU+  DTH**2) 

TERM2*(2.*  ATH(J))  /  <<<  l.-ATH(J)**2  )**2  )*  FSDTH) 

ALAM ( J ) =  TERM1  *  TERM2 
400  CONTINUE 
JLONG=360 

-  COMPUTE  CTH  =  F ( CLAM)  *  WHERE 

-  CLAM  =  10**L  »  BY  USING  ALAM*  F(ATH)  TABLE 

425  DO  750  J=I0FL1*I0FL2 *IOFL3 
CL ( J  ) *  BL( J ) 

CLAMU)*  BLAM(J) 

L ( 1 ) *  XLOCF I  L ( 1 ) ) 

L ( 2 ) *  XLOCF ( ALAM(1 ) ) 
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L(3)=  XLOCF ( ATH ( 1 )  ) 

L<4>=1 
L  ( 5 )  =  1 
L ( 6 ) =5 
L  ( 7 )  =  JLONG 
L ( 8 ) =0 

CTH( J) =  TABF(CLAM( J| *L( 1) I 
GO  TO  (750,710  ),I 
710  PRINT  720,  IOFD,  J,CTH(J) 

720  FORMAT (16H1CTH=F( CLAM)  BAD2 1 5 *E20«8 ) 

GO  TO  9999 

750  CONTINUE 

-  COMPUTE  CRHOT 

DO  800  J=I0FL1 » IOFL2, IOFL3 
CRHOT ( J ) =  CTH(J)**2  /  (1.-  CTH(J)**2  ) 
800  CONTINUE 

-  COMPUTE  CRHO  BY  NEWTON  RAPHSON 

DO  900  J=IOFLl,IOFL2*IOFL3 
RHONEW*CRHOT ( J ) 

AGPRN*0. 

C 

DO  880  IT=1,MAXIT 
I F  C IT-10)825, 825*830 
825  RHO*RHONEW 
GO  TO  850 
830  RHO=RHONEW 

IF (AGPRO-AGPRN) 835,835 ,840 
835  GO  TO  890 
840  GO  TO  850 
C 

850  AGPRO=AGPRN 
V*D( IOFD) 

FLAM*10.**CL< J) 

THETA*SQRTF(RH0/(1.+RH0) ) 

CON* 1 . / ( SQRTF ( FLNU ) *BETA ) 
EXP=(-,5)*(FLNU+1. > 

DTH=V/THETA 

FSDTH=CON  *  ( l.+DTH**2/FLNU)**EXP 
Tl*( *5*FLAM)  /  ( SQRTF (RHO*(RHO+l. )**3) ) 
T2=(FLNU+DTH**2)  /  (FLNU-1.) 

GPR* ( T 1*T2*F5DTH) -1 # 

T A=FLAM/  ( 4** ( FLNU-1 • ) ) 

TB=SQRTF ( RHO* ( RHO+1 . ) ) 
TQA«-4.*(FLNU+V**2)*RHO**2 
TQB«(V**2*FLNU-6.*V**2-FLNU)  *  RHO 
TQC=V**2*(FLNU-2.) 

U1*TA/TB 

U2= ( TQA+TQB+TQC ) / ( TB**2 ) 
G2PR*(Ul*U2*FSDTH) /(TB**2) 

RHONEW*RHO  -  ( GPR/G2PR ) 

AGPRN*ABSF ( GPR ) 


“  32  - 


Program  D.1,  page  7 


880 

C 

890 


C 

900 

C 

c 


910 

920 

950 

C 

C - 

C 


C 

D960 

D 

D 

C 

C 

D 

D 

D 

C 

C - 

C 

D 

D 

D 

D 

C 

D 

D 

D 

D 

D 

D970 

975 

C 

D980 

D 

C 

C 

D 

C 

D 


CONTINUE 

CRHO( J) “RHONEW 
ITER ( J ) =IT 
CGPR ( J ) XGPR 

CONTINUE 

COMPUTE  CLGRHO 

DO  950  J»IOFLl*IOF.L2*IOFL3 
IF(CRH0(J)-0. 1910 *910 *920 
CLGRHO! IOFD*J) =-99*9999 
GO  TO  950 

CLGRHOI IOFD* J)=L0G10F (CRHOI J) ) 

CONTINUE 

COMPUTE  CNG*  CLGNG 

DO  1100  J=I0FL1 » IOFL2 » I0FL3 
DTH»D( IOFD) /CTH ( J ) 

IF (DTH-8.) 1010*960*960 

TNU*FLNU 

TBETA=BETA 

H*DTH 

DENSITY  FST 

2= (l./SQRTF(TNU ))*(!. /BETA! 

EXP=T-»5 )*( TNU  +  1*  > 

FSTxZ*(l.+H**2/TNU)»*EXP 

RIGHT  TAIL  GST 

CONA-TNU/ ( TNU+T**2 ) 

TERM=1./TNU 
SUM=  l./TNU 
SUMB4=0  « 

DO  975  JK*2  *500 
FJ=JK 

TINCR=((2.*FJ-3.)*(TNU+2**FJ-4.) )/( <2.*FJ-2* >*<TNU+2.*FJ-2. ) ) 

TERM=TINCR*CONA*TERM 

SUM“SUM+TERM 

IF (SUM-SUMB4) 980*980*970 

SUMB4*SUM 

CONTINUE 

CONB=< l./TBETA)*CONA**(.5*TNU) 

GST=CONB*SUM 

LINEAR  LOSS  FLST 

FLST*( (TNU+H**2)/(TNU-1. >)*FST-H*GST 
FLSDTH»FLST 
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GO  TO  1050 
C 

1010  L ( 1 ) *  XLOCF (  L ( 1 )  ) 

L ( 2  >  *  XLOCF  ( T ( 1 )  ) 

L  ( 3 )  =  XLOCF ( FLS ( 1 )  ) 

L  ( 4 )  *  1 
L  (  5  )  =  1 
L(6)  =  5 
L  (  7 )  =  1001 
L(8)»  0 

FLSDTH  =  T ABF (  DTH*L( 1 ) ) 

I=L( 8) 

GO  TO  ( 1050*1020) ♦  I 
1020  PRINT  1030*  IOFD*  J»  FLSDTH 
1030  FORMAT ( 15H1FLSDTH  BAD2 1 5 *E20.8 ) 

GO  TO  9999 
C 

1050  CLAMTT *CLAM( J ) 

D  CLAMT*CLAMTT 

CTHTT*CTH(J) 

D  CTHT  *CTHTT 

CRHOTT  »CRHO ( J ) 

D  CRHOT*CRHOTT 

D  CNG( J)»CLAMT*CTHT*FLSDTH-CRHOT 

IF(CNG(J)-0. ) 1060*1060*1070 
1060  CLGNG ( IOFD, J ) -  -99.9999 
GO  TO  1100 

1070  CLGNG ( IOFD, J ) *  LOGIOF  (CNG(J)) 

C 

1100  CONTINUE 
C 

C  -  RHO  ASYMPTOTIC  CALCULATION  FOR  VARIOUS  CLAM  VALUES 

C 

KOFD=(10*I0FD)-9 

FSOFD=FS(KOFD) 

C 

CONST* { .5*( ( FLNU+DI IOFD) #*2 ) / (FLNU-1 • ) )*FSOFD )**.5 
C 

DO  1600  J=I0FL1 * IOFL2  » IOFL3 
1600  CRHOASt  J ) =CONST*CLAM( J)**.5 
C 

DO  1650  J*I0FL1 » IOFL2  * IOFL3 
1650  RHORAT ( J )=CRHO( J ) /CRHOAS ( J ) 

C 

C  -  WRITE  MAIN  FUNCTIONS  ON  TAPE  6 

C 

WRITE  OUTPUT  TAPE  6»81*D(IOFD) 

WRITE  OUTPUT  TAPE  6*61 
WRITE  OUTPUT  TAPE  6*22 

WRITE  OUTPUT  TAPE  6*71  * ( CL ( J ) *CLAM( J ) ,CTH( J ) ,CRHO( J ) *CLGRHO( IOFD  * 
XJ),CNG< J)»CLGNG( IOFD* J) , J* IOFL1 * IOFL2* IOFL3 ) 

C 

C  -  WRITE  CLGRHO  AND  CLGNG  ON  TAPE  2 

C 

WRITE  OUTPUT  TAPE  2 » 15 * ( CLGRHO ( IOFD* J) • J* IOFL1 » IOFL2  * IOFL3 ) 

WRITE  OUTPUT  TAPE  2 » 15* ( CLGNG! IOFD, J ) * J-IOFL1 , IOFL2 » IOFL3 ) 

C 
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-  WRITE  CRHO  FUNCTIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6»81*D(IOFO) 

WRITE  OUTPUT  TAPE  6»22 

WRITE  OUTPUT  TAPE  6 .72  * ( CL ( J ) »CRHOT ( J ) »CRHO< J ) »CGPR ( J ) • I  TER (J  )  * J= I 
X0FL1 *I0FL2»I0FL3) 

-  WRITE  ASYMPTOTIC  CALCULATIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6»81»D(IOFD) 

WRITE  OUTPUT  TAPE  6.63 
WRITE  OUTPUT  TAPE  6*22 

WRITE  OUTPUT  TAPE  6 »73* ( CL ( J ) »CRHO( J ) tCRHOAS (J ) tRHORAT ( J ) » J= IOFL 1 
X  » IOFL2  » IOFL3 ) 

-  COMPUTE  PHI  AND  PSI  FOR  D 

PHIASI IOFD)*. 5  *  LOG 1 OF (  .5*  ( (FLNU+DI IOFD) **2 ) / ( FLNU-1 • ) ) *FSOFD) 

FLSOFD=FLS<KOFD> 

PS  IASI IOFD) = LOG 1 OF (FLSOFD) 

PSIASl IOFD)=LOG10F (FLSOFD) 

-  INDEX  ON  D 

2000  CONTINUE 

-  WRITE  SS  AND  NG*  PLUS  PHIAS  AND  PSIAS 

CALL. WRSSNG ( NU*CL  »CLGRHO»CLGNG»PH IAS»PSIAS) 

-  INDEX  ON  NU 

3000  CONTINUE 
C 

9999  CONTINUE 

END  FILE  2 
END  FILE  8 
REWIND  2 
REWIND  8 
PRINT  6 
PAUSE 
CALL  EXIT 
END 
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D.2  R  and  P  for  3  <  V  <  00 


-  OPTIMAL  SAMPLE  SIZE  AND  NET  GAIN 


PRINT  5 
PAUSE 

DIMENSION  CL ( 121 ) *CLGRHO (41*121) »CLGNG( 41*121  ) 
COMMON  CL*CLGRHO*CLGNG 
DIMENSION  PH I  AS ( 12 1 ) »PS  I AS (121) 

COMMON  PH I  AS. PS  I AS 


C 

5 

6 

10 

15 

16 
21 
22 

51 

52 

53 
61 

63 

71 

72 

73 
81 

C 

C 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 


DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 

DIMENSION 


FLNUT ( 19 ) »STVT ( 19 ) .BETAT ( 19 ) 

T ( 1001) .  FS( 1001 ) ,GS( 1001 ) »FLS< 1001 ) 

T’S  ( 91 ) »  FSS ( 91  ) 

L(  11 ) 

D(  41 ) .THCTT (  41 ) »THCT (  41)»FLAMCT(  41)»FLCT(  41) 
ATH(360>  *BL(121> .BLAM(121) 

ALAM < 360 ) .CLAM ( 121) »CTH< 121) 

CRHOT  (121)  *CRHO(  121). ITER  (121 )  .CGPRU21 ) 

CRHOAS( 121 ) .RHORAT (121) 

CNGI121) 


FORMAT ( 38H1MOUNT  OUTPUT  TAPES  ON  LOGICAL  2  AND  8) 
FORMAT ( 25H1DI SMOUNT  LOGICAL  2  AND  8) 

FORMAT (31 10) 

FORMAT ( 5E20.8 ) 

FORMAT ( 5F20»8 ) 

FORMAT ( 1H1 ) 

FORMAT ( 1H  ) 

FORMAT (21H0STUDENT 
FORMAT ( 24H0STUDENT 
FORMAT ( 25H0STUDENT 
FORMAT ( 110H0 

L  CRHO  CLGRHO  CNG 


FORMAT (70H0 


DENSITY.  NU= 14 ) 

RIGHT  TAIL.  NU=I4) 
LINEAR  LOSS.  NU=I4) 

CL  CLAM 

CLGRHO 

CL  CRHO 


CTH 

CLGNG) 

CRHOAS 


RHORAT ) 


FORMAT(F10.2*3E20.8.F10.4,E20.8»F10.4) 
FORMAT (F10*2.3E20«8.I10) 

FORMAT (F10*2.3E20*8) 

FORMAT ( 3H0D=F6»2 ) 


SPECIFY  FLNUT.  STVT.  AND  BETAT 


FLNUT ( 4 ) =3» 
FLNUT ( 5 ) =4« 
FLNUT ( 6 ) =5 • 
FLNUT ( 7 ) -6« 
FLNUT ( 8 ) =7 • 
FLNUT ( 9 ) -8 • 
FLNUT ( 1 0 ) =9 • 
FLNUT ( 1 1 ) =10« 
FLNUT ( 12 ) =12  • 
FLNUT ( 13 ) =15» 
FLNUT ( 14) =20. 
FLNUT { 15 ) =24» 
FLNUT ( 16 ) =30 • 
FLNUT ( 17 ) =40 • 
FLNUT ( 18 ) =60 • 
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FLNUT( 19)=120* 

STVT ( 4 )  =  0.1064199529207075E-02 
STVT ( 5 ) =  0* 281001 81 13579956 E-03 
STVT (6 ) =  0.8547378787148180E-04 
STVT (7 )  =  0* 289 599 1377476 8 13 E-04 
STVT ( 8 ) =  0* 10697 10 1445 3864 IE— 04 
STVT  ( 9 )  =  0  •  42440907638 14246 E**0 5 
STVT ( 10 ) =0 • 1789 11 87 15962 368 E-05 
STVT ( 1 1 ) =0 . 7947765877982060E-06 
STVT (12) =0. 1790661843836 187E-06 
STVT ( 13 ) =0 » 2498 449 07 146 2472 E— 07 
STVT (14) =0.1 581890879357 194E-08 
STVT ( 15 ) =0 • 2 45 7 762 76 383 7653 E— 09 
STVT ( 16 ) =0*22876257041 14807E- 10 
STVT (17) =0 .9656558502057788 E— 12 
STVT ( 18 ) =0 . 1068842957337041 E-13 
STVT ( 19 ) =0  *856972 05 876 805 08 E— 17 

BETAT (4 ) =  0.1570796326794897E+01 
BETAT ( 5 )  =  0.1333333333333333E+01 
BETAT ( 6 )  =  0. 11 7809 7 24 5096 17 2 E+01 
BETAT ( 7 ) =  0 • 106666666666666 7 E+01 
BETAT ( 8 )  =  0 • 98 1 747 704246 81 04E+00 
BETAT ( 9 )  =  0.9142857142 857 143 E+00 
BETAT ( 10 ) =0.8 590 29241 2 159591 E+00 
BETAT( 11 ) =0*8 126984 126984 12 7 E+00 
BETAT ( 12 )=0*7388167388167388E+00 
BETAT ( 13 >=0.6580777580029401 E+00 
BETAT ( 14 )=0.5675463855030419E+00 
BETAT ( 15) =0.5 170194816176779E+00 
BETAT ( 16) =0. 461474553400974 1E+00 
BETAT ( 17 )=0 *39881 730689488 10E+00 
BETAT ( 18 ) =0 • 3249554203948302E+00 
BETAT ( 19 )=0.2293000137934539E+00 

-  GENERATE  T  AND  TS 

I  T  =  0 
T ( 1 ) =IT 

DO  110  1=2*1001 
I T=  IT  +  1 
TT=IT 

110  T ( I )  =  TT/100. 

IT  =  9910 
TS( 1 )=  IT/1000 
DO  115  1=2*91 

IT=IT+1 
TT=IT 

115  TS( I )=TT/1000. 

- D  LIST  0  <  .1)4 

ID  =  0 
D ( 1 ) =  ID 
DO  130  1=2.41 
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ID=ID+1 
TD=  ID 

130  D ( I ) *  TD/10. 

-  ATH  LIST  360  VALUES  .9999  9999  TO  .01 


ATH( 

1 ) *.9999 

9999 

ATHI 

2)*=. 9999 

9998 

ATH  ( 

3 ) =.9999 

9997 

ATH( 

4 ) =.9999 

9996 

ATHI 

5 ) =.9999 

9995 

ATHI 

6 ) =.9999 

9994 

ATH  ( 

7 ) =.9999 

9993 

ATHI 

8 ) =.9999 

9992 

ATHI 

9 )=.9999 

9991 

ATHI 

10)=. 9999 

999 

ATHI 

11)=. 9999 

998 

ATHI 

12)=. 9999 

997 

ATHI 

13)=. 9999 

996 

ATHI 

14)=. 9999 

995 

ATHI 

15)=. 9999 

994 

ATHI 

16)=. 9999 

993 

ATHI 

17)=. 9999 

992 

ATHI 

18)=. 9999 

991 

ATHI 

19)=. 9999 

99 

ATHI 

20)=. 9999 

98 

ATHI 

21)=. 9999 

97 

ATHI 

22)=. 9999 

96 

ATHI 

23)=. 9999 

95 

ATHI 

24)=. 9999 

94 

ATHI' 

25)=. 9999 

93 

ATHI 

26)=. 9999 

92 

ATHI 

27)=. 9999 

91 

ATHI 

28)=. 9999 

9 

ATHI 

29)=. 9999 

8 

ATHI 

30)=. 9999 

7 

ATHI 

31 )  =.9999 

6 

ATHI 

32) =.9999 

5 

ATHI 

33  >  =  .9999 

4 

ATHI 

34)=. 9999 

3 

ATHI 

35)=. 9999 

2 

ATHI 

36)=. 9999 

1 

ATHI 

37)=. 9999 

DO  133  J=38 * 236 
133  ATHI J)=ATHl J-l )-. 0001 
ATH( 237 ) =.979 
DO  134  J=238.266 

134  ATH(  J)=ATH(  J-D-.001 
ATH( 267) =.94 

DO  135  J«268»360 

135  ATH< J)=ATH( J-l ) — *0 1 

WRITE  OUTPUT  TAPE  6 . 15  *  I ATH I J )  ♦  J= 1 . 360 ) 
WRITE  OUTPUT  TAPE  6  *  16 . I ATH I J ) ♦ J= 1 . 360  ) 

-  BL  AND  BLAM  L=-2.(.l)10 
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LT  *-20 
BL(l>»-2. 

BLAM ( 1 ) =*01 
DO  160  J=2* 121 
LT*LT+1 

1F(LT  -  0)155*150.155 
150  BL(J)»  +  0 
GO  TO  160 
155  FLT  *  LT 

BL ( J  )*  FLT/10. 

160  BLAM(J)  *  10.**  BL(J) 

—  READ  NUMNU  »  MAIN  INDEX  ON  NU  -  IOFNU 

READ  INPUT  TAPE  5,10*NUMNU 
WRITE  OUTPUT  TAPE  6»10»NUMNU 

DO  3000  IOFNU  =  1 » NUMNU 

-  READ  INU.  NU*  AND  MAXIT 

READ  INPUT  TAPE  5*10.INU*NU*MAXIT 
WRITE  OUTPUT  TAPE  6»10*INU»NU»MAXIT 

FLNU=FLNUT (INU) 

STV*STVT( INU) 

BETA*BETAT( INU) 

WRITE  OUTPUT  TAPE  6.15»FLNU.STV.BETA 

-  INITIALIZE  CLGRHO  AND  CLGNG 

DO  200  1=1*41 
DO  190  J=1 *121 

CLGRHO(I.J)  =  0. 

CLGNG ( I  * J)  *0. 

190  CONTINUE 
200  CONTINUE 

-  COMPUTE  FS 

TERM1=  1 •/( BET A*SQRTF (FLNU )  ) 

EXP=-( (FLNU+1. )/2. ) 

DO  230  1=1*1001 

230  FS  ( I ) =TERM1*  (<1.+(T  (I)**2 ) /FLNU)**EXP) 

-  COMPUTE  GS 

GS(lOOl)  =  STV 
DO  235  1=1*91 

235  FSS( I ) =TERMl*  ( 1 1 .  +  ( TS ( I ) **2 ) /FLNU) **EXP ) 

DO  250  1=2*10 
J=  11-1 
X«  10*(J-1) 

D  Tl=16067  *  * ( FSS( X+l )+FSS( K+ll ) ) 
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T2*106300.  * ( FSS( K+2 )+FSS( K+10) ) 
T3*48525*  * ( FSS ( K+3 )+FSS< K+9 ) > 

T 4*272400*  *<FSS(K+4)+F$S(K+8 |) 
T5-260550.  * ( FSS(IC+5  )+FSS(  K+7  » ) 

T 6*427368 .  * ( FSS ( K+6 >  ) 

ATI* C 5. *0.001 1/299376. 
AT2*T1+T2-T3+T4-T5+T6 
A*AT1*AT2 
L*1002-I 
M*1002-I+l 
230  GS(L)«GS<M)+A 

DO  260  1*11^1001 
J«  1001-1 

Tl«16067.  *( FS  (J+1I+FS  <J+11>) 

T 2* 106300*  »(FS  (J+2H-FS  (J+10)) 
T3-48525.  *(FS  (J+3I+FS  <J+9>> 
T4*272400.  *(FS  (J+4KFS  <J+8>> 
T5-260550.  *(FS  (J+5I+FS  (J+7)) 

T 6*427368 •  *(FS  (J+6)) 
AT1*<5.*0.01  )/299376. 
AT2*T1+T2-T3+T4-T5+T6 
A*AT1*AT2 
K»1002-I 
L»1002-I+10 
260  GS(K>*GS(LM-A 

-  COMPUTE  FLS 


BOT  =  FLNU-1. 

DO  270  1*1*1001 

270  FLS{ I )* ( (FLNU+  T ( I >**2. ) /BOT )*FS( I )  -  T  ( I  )* 

-  PRINT  STUDENT  DISTRIBUTIONS 

WRITE  OUTPUT  TAPE  6*51 *NU 

WRITE  OUTPUT  TAPE  6*15»<FS  (I) » 1*1 . 1001 *100 ) 
WRITE  OUTPUT  TAPE  6*52*NU 

WRITE  OUTPUT  TAPE  6»15*(GS  (I ) *  1*1  *  1001 *100 ) 
WRITE  OUTPUT  TAPE  6*53*NU 

WRITE  OUTPUT  TAPE  6*15* ( FLS( I) *  I *1* 1001 *100 ) 

-  READ  INDEX  IOFD  INDEX  ON  D 

READ  INPUT  TAPE  5* 10* I0FD1 * I0FD2. I0FD3 
WRITE  OUTPUT  TAPE  6.10.I0FD1.1OFD2.IOFD? 

DO  2000  IOFD  *  I0FD1 * I0FD2  » I0FD3 

-  READ  AND  PRINT  THCT ( IOFD) *  IOFLl . I0FL2. I0FL3 

READ  INPUT  TAPE  5* 15 *THCT (IOFD ) 

READ  INPUT  TAPE  5 *  10* IOFLl * I0FL2* I0FL3 

WRITE  OUTPUT  TAPE  6*21 

WRITE  OUTPUT  TAPE  6»15.THCT< IOFD> 

WRITE  OUTPUT  TAPE  6.10. IOFLl ♦ I0FL2 . I0FL3 


GS<  I) 
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C 

C  -  COMPUTE  ALAM  *  FCATH)  WHERE  ATH  *  OR  GREATER  THAN  THCT ( IOFD) 

C  -  NOTE  ATH  IS  STORED  IN  DECREASING  ORDER 

C 

DO  400  J*1.360 

IF ( ATH( J I-  THCTC IOFD) >350,350.370 
C 

350  JLONG*J— 1 
GO  TO  425 
C 

370  DTH  *  DC IOFD) /  ATH(J) 

LCD*  XLOCFCLC1) >  ' 

LC2)*  XLOCFCTCl)) 

L ( 3 ) *  XLOCF ( FSC 1 ) ) 

L(4) *  1 
L  ( 5 )  *  1 
L  ( 6 )  *  5 
L ( 7 ) *  1001 
L  ( 8 )  *  0 

FSDTH  *  TABFCDTH.LC1) > 

I *LC  8 ) 

GO  TOC  390*375  )  *  I 
C 

375  PRINT  380*  IOFD*  J»  DTH 

380  FORMATC 15H1TFSC D/ATH)  BAD2I5*E20.8) 

GO  TO  9999 
C 

390  TERM1*C  FLNU-1* ) /CFLNU+  DTH**2) 

TERM2*( 2«*  ATHC J) )  /  CCC  i.-ATHCJ)**2  )**2  )*  FSDTH) 

ALAMC J) *  TERM1  *  TERM2 
400  CONTINUE 
JLONG*360 
C 

C  -  COMPUTE  CTH  *  FCCLAM)  *  WHERE 

C  -  CLAM  *  10**L  ♦  BY  USING  ALAM*  FCATH)  TABLE 

C 

425  DO  750  J»IOFLl*IOFL2*IOFL3 
CL C J ) *  BLCJ) 

CLAM  C J ) *  BLAMC J ) 

LCD*  XLOCF C  LCD) 

L C 2 ) *  XLOCF « ALAMC 1 ) ) 

L ( 3 ) *  XLOCF  C  ATH ( 1 ) ) 

L  <  4 ) *1 
L  ( 5 )  *1 
L(6)  *5 
L  <  7 ) *  JLONG 
L  C  8 )  *0 

CTHCJ)*  TABFCCLAMC J) *L( 1) ) 

GO  TO  C 750*710  )*I 
710  PRINT  720.  IOFD.  J. CTHCJ) 

720  FORMATC 16H1CTH*FC CLAM)  BAD2I5.E20.8 ) 

GO  TO  9999 
C 

750  CONTINUE 
C 

C  —  COMPUTE  CRHOT 
C 
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DO  800  J«10FL1?>I0FL2»I0FL3 
CRHOT ( J) «  CTH(J)**2  /  <1.-  CTH(J)**2  V 
800  CONTINUE 

—  COMPUTE  CRHO  BY  NEWTON  RAPHSON 

DO  900  J*I0FL1.I0FL2» IOFL3 
RH0NEW*CRH0T< J) 

AGPRN*0. 

DO  880  IT*1 .MAXIT 
I F  C IT— 10)825.825*830 
825  RHO-RHONEW 
GO  TO  850 
830  RHO*RHONEW 

I F ( AGPRO-AGPRN ) 835 *835 »840 
835  GO  TO  890 
840  GO  TO  850 

850  AGPRO*AGPRN 
V*DUOFD> 

FLAM*10.**CL(J) 

T  HET  A*SQRT  F ( RHO/ ( 1 • +RHO ) ) 

CON« 1 . / ( SORT F ( FLNU ) *BET A I 
EXP*  C— #5 )*  <  FLNU+1 • ) 

DTH*V/THETA 

FSDTH-CON  *  <1.+DTH**2/FLNU)**EXP 
Tl=*  ( •5*FLAM)  /  (5QRTF (RHO* ( RHO* 1* )**3) ) 
T2*(FLNU+DTH**2>  /  (FLNU-1.) 
GPR»(T1*T2*FSDTH)-1. 

TA-FLAM/  (4.*<FLNU-1.>) 

TB-SQRTF ( RHO* ( RHO* 1 . )  ) 

TQA*-4. * ( FLNU+V**2 ) *RH0**2 
TQB»(V**2*FLNU-6.*V**2-FLNU)  *  RHO 
TQC«V**2*(FLNU-2.) 

U1*TA/TB 

U2» ( TQA+TQB+TQC > / ( TB**2 ) 
G2PR»(U1*U2*FSDTH)/(TB**2) 

RHONEW*RHO  -  (GPR/G2PR) 

AGPRN*ABSF ( GPR ) 

880  CONTINUE 

890  CRHO ( J ) *RHONEW 
ITER ( J ) *1 T 
CGPR(J)*GPR 

0  CONTINUE 

—  COMPUTE  CLGRHO 

DO  950  J*I0FL1.I0FL2.I0FL3 
IF(CRHO(J)-0. >910.910.920 
910  CLGRHOI IOFD»J)«- 99.9999 
GO  TO  950 

920  CLGRHOI IOFD. J ) *LOG10F ICRHOl J )) 

950  CONTINUE 


-  42  - 


nnn  nnn  nnn  n  n  n  nnn  n  n  non 


Program  D.2,  page  8 


-  COMPUTE  CNG»  CLGNG 

1000  DO  1100  J*  I0FL1 * I0FL2 * I0FL3 
DTH*  D(IOFD)/  CTH(JJ 
L ( 1 ) *  XL0CF (  LCD) 

L  ( 2 )  *  XLOCF  ( T .( 1 ) ) 

L<3)«  XLOCF (FLSC 1 )  ) 

LCD*  1 
L  C  5 )  *  1 
L  C  6 )  *  5 
L  C  7) *  1001 
LC8)*  0 

FLSDTH  *  TABFC  DTH*LC 1 ) ) 
i*L<e> 

GO  TO  (1050*1020).  I 
1020  PRINT  1030*  IOFO*  J.  FLSDTH 
1030  FORMAT (15H1FLSDTH  BAD215,E20.8) 

GO  TO  9999 

1050  CNGCJ)*  CLAM(J)  *  CTHIJJ  *  FLSOTH  -  CRHOCJ) 

IFCCNGCJ1-0*) 1060* 1060*1070 
1060  CLGNG( I0FD* J )*  -99.9999 
GO  TO  1100 

1070  CLGNG ( IOFD*J)«  LOGIOF  CCNGCJ)) 

1100  CONTINUE 

-  RHO  ASYMPTOTIC  CALCULATION  FOR  VARIOUS  CLAM  VALUES 

KOFD*(10*IOFD>-9 

FSOFD«FS(KOFD) 

CONST- C . 5*  C ( FLNU+D ( I OFD ) **2 ) / C  FLNU-1 • > ) *FSOFD )  **.  5 

DO  1600  J*I0FL1 * IOFL2 * IOFL3 
1600  CRHOAS  C  J ) “CONST *CLAM ( J  >  ##• 5 

DO  1650  J-IOFLl » IOFL2  * IOFL3 
1650  RHORATC J) -CRHOCJ J/CRHOASCJ) 

-  WRITE  MAIN  FUNCTIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6*B1»DCI0FD) 

WRITE  OUTPUT  TAPE  6.61 
WRITE  OUTPUT  TAPE  6.22 

WRITE  OUTPUT  TAPE  6.71. (CLCJ) .CLAMCJ).CTH(J) .CRHOCJ) .CLGRHOC IOFD 
XJ). CNGCJ) . CLGNG ( IOFD. J ) ♦ J* IOFL1 ♦ IOFL2. 10FL3 ) 

-  WRITE  CLGRHO  AND  CLGNG  ON  TAPE  2 

WRITE  OUTPUT  TAPE  2 .15. < CLGRHO ( IOFD. J ) .J-IOFL1 » IOFL2 . IOFL3 ) 

WRITE  OUTPUT  TAPE  2 .15* C CLGNG l IOFD. J) . J- IOFL1 ♦ IOFL2. IOFL3 ) 

-  WRITE  CRHO  FUNCTIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6. 81. DC IOFD) 
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WRITE  OUTPUT  TAPE  6*22 

WRITE  OUTPUT  TAPE  6»72 » ( CL< J) tCRHOT ( J) »CRHO<J J »CGPR( J) *  ITER ( J ) ♦ J« I 
X0FL1 » IOFL2  * I0FL3 ) 

-  WRITE  ASYMPTOTIC  CALCULATIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6*81»D(I0FD) 

WRITE  OUTPUT  TAPE  6»63 
WRITE  OUTPUT  TAPE  6.22 

WRITE  OUTPUT  TAPE  6»73, ( CL (J I »CRHO< J) .CRHOAS(J) .RHORAT ( J ) * J» IOFL1 
X* I0FL2 ♦ I0FL3 ) 

—  COMPUTE  PHI  AND  PS  I  FOR  D 

PHIASI IOFD)*.5  *  LOG 1 OF ( «5* ( (FLNU+OI IOFD)**2 ) / (FLNU”1» ) I*FSOFD) 

FLSOFD*FLS  ( ICOFD ) 

PSIASI IOFD)»LOG10F (FLSOFDI 
PSIASl IOFD)*LOG10F (FLSOFO) 

-  INDEX  ON  D 

2000  CONTINUE 

-  WRITE  SS  AND  NG*  PLUS  PHIAS  AND  PSIAS 

CALL  WRSSNG ( NU • CL » C  LGRHO »  CLG  NG  »  PH I AS  *  PS I AS ) 

-  INDEX  ON  NU 

3000  CONTINUE 
C 

9999  CONTINUE 

END  FILE  2 
END  FILE  8 
REWIND  2 
REWIND  8 
PRINT  6 
PAUSE 
CALL  EXIT 
END 
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C  -  OPTIMAL  SAMPLE  SIZE  ANO  NET  GAIN*  NU«INF  '• 

c 

PRINT  5 
PAUSE 
C 

DIMENSION  CL ( 121 ) *CLGRHO( 41 *121 ) *CLGNG(41 *121 ) 

COMMON  CL*CLGRHO*CLGNG 
DIMENSION  PHI  AS ( 121 ) *PSI AS! 121 ) 

COMMON  PHIAS.PSIAS 
C 

D  DIMENSION  T< 1001)*  FS< 1Q01 ) *GS< 1Q01 ) *FL$< 1001 ) 

D  DIMENSION  TSI91 ) *FSS<91 ) 

DIMENSION  L(ll) 

DIMENSION  D<  41)»THCTT!  41)*THCT<  41)*FLAMCT!  41 ) »FLCT (  41) 
DIMENSION  ATH( 360) *BLI121) *BLAM( 121) . 

DIMENSION  ALAM( 360 ) *CLAM ! 121 ) *CTH! 121 ) 

DIMENSION  CRHOT 1121) .CRHO!  121 )  .  I TERU21 )  ,CGPR  <121  > 

DIMENSION  CNG ( 121 ) 

DIMENSION  CRHOASI 121) tRHORAT (121) 

C 

5  FORMAT (38H1M0UNT  OUTPUT  TAPES  ON  LOGICAL  2  AND  8) 

6  FORMAT ( 25H1DISMOUNT  LOGICAL  2  AND  8) 

10  FORMAT (3110) 

15  FORMAT! 5E20. 8) 

16  FORMAT! 5F20, 8) 

21  FORMAT (1H1) 

22  FORMAT! 1H  ) 

51  FORMAT! 2 1HOSTUDENT  DENSITY*  NU*I4) 

52  FORMAT t 24H0STUDENT  RIGHT  TAIL*  NU«14) 

53  FORMAT! 25H0STUDENT  LINEAR  LOSS*  NU«I4) 

61  FORMAT  U10H0  CL  CLAM  CTH 

1  CRHO  CLGRHO  CNG  CLGNG) 

63  FORMAT! 70H0  CL  '  CRHO  CRHOAS 

X  RHORAT) 

71  FORMAT (F10*2*3E20*8*F10«4*E20*8*Fl6«4) 

72  FORMAT!F10*2»3E20. 8*110) 

73  FORMAT (F10V2*3E20*8) 

81  FORMAT! 3H0D=F6#2) 

C  ' 

C  -  NU=999 

C  • 

NU»999 

C 

C  —  GENERATE  T  AND  TS 
C 

IT*0 

D  T1 1 )*1T 

DO  110  1*2*1001 
IT*  IT  +  1 

D  TT*IT 

D  110  T ! I )  *  TT/100. 

C 

IT  *  9910 

D  TS(1 )*  IT/1000 

DO  115  1=2*91 

IT*IT+1 
D  TT*I T 
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11  5  T${ 11*77/1000. 

— -  D  LIST  oum 

ID  *  0 
D(l)=  ID 
DO  130  1*2,41 
ID*ID+1 
TD=  ID 

130  DU1*  TD/10. 

-  ATH  LIST  360  VALUES  .9999  9999  TO  .01 


ATHC 

1 )*»9999 

9999 

ATH( 

2 >*.9999 

9998 

ATH( 

3)*. 9999 

9997 

ATH( 

4)*. 9999 

9996 

ATH( 

5 >=.9999 

9995 

ATHJ 

61**9999 

9994 

ATH( 

71**9999 

9993 

ATH< 

8 1 **9999 

9992 

ATHJ 

9)  *,9999 

9991 

ATHJ 

101**9999 

999 

ATHJ 

111**9999 

998 

ATHJ 

121**9999 

997 

ATHJ 

131**9999 

996 

ATHJ 

141**9999 

995 

ATHJ 

151**9999 

994 

ATHJ 

161**9999 

993 

ATHJ 

171**9999 

992 

ATHJ' 

181**9999 

991 

ATHJ 

191**9999 

99 

ATHJ 

201**9999 

98 

ATHJ 

211**9999 

97 

ATHJ 

221**9999 

96 

ATHJ 

231**9999 

95 

ATHJ 

241**9999 

94 

ATHJ 

251**9999 

93 

ATHJ 

261**9999 

92 

ATHJ 

271**9999 

91 

ATHJ 

281**9999 

9 

ATHJ 

291**9999 

8 

ATHJ 

301**9999 

7 

ATHJ 

311**9999 

6 

ATHJ 

3.31**9999 

5 

ATHJ 

331**9999 

4 

ATHJ 

341**9999 

3 

ATHJ 

351**9999 

2 

ATHJ 

361**9999 

1 

ATHJ 

371**9999 

DO  133  J*38»236 

133  ATHJ J)*ATHJJ~1>~. 0001 
ATM  237 )  *.979 

DO  134  J“238  *266 

134  ATHJ  J)»ATHI  J-11-.001 
ATHJ 267 1 *»94 

DO  135  J*268  *360 
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135  ATH(  J)*ATH(  J-D-,01 

c 

WRITE  OUTPUT  TAPE  6*15* ( ATH( J) »J* 1*360) 
WRITE  OUTPUT  TAPE  6* 16 » ( ATH ( J ) * J*1 »360 ) 
C 

C  -  BL  AND  BLAH  L*-2*(.l)10 

C 

LT*-20 

BL(l)*-2* 

BLAM ( 1 ) **01 
DO  160  J*2»121 
LT*LT+1 

IF(LT  -  0) 155*150*155 
150  BL(J)«  +  0 
GO  TO  160 
155  FLT  *  LT 

BL( J  )*  FLT/10. 

160  BLAM(J)  *  10***  BL(J) 

C 

C  -  INITIALIZE  CLGRHO  AND  CLGNG 

C 

DO  200  1*1 *41 

DO  190  J*l*121 

CLGRHOI I *J)  *  0* 

CLGNG ( I  * J )  *0* 

190  CONTINUE 
200  CONTINUE 
C 

C  -  COMPUTE  FS 

C 

D  CON*l./SQRTF(2.*3. 141592653589793) 

C 

DO  230  1*1*1001 

D  230  FS< I )*CON*EXPF( ( — • 5 1 *T ( I »  I 
C 

C  -  COMPUTE  GS 

C 

D  GStlOOl)*. 76198  53024  16059  16  E-23 

C 

DO  235  1*1*91 

D  235  FSS( I )*CON  *  EXPFI (-.5)*TS(I )»*2) 

C 

DO  250  1*2.10 
J*  11-1 
K*  10*(J-1) 

D  Tl*16067*  * ( FSSI K+l )+FSS( K+ll ) ) 

D  T2*106300.  *(FSS(K+2)+FSS(K+10) ) 

D  T3*48525.  *(FSSU+3)+FSS(  K+9)  ) 

D  T4*272400.  *(FSS(K+4)+FSS< K+8) ) 

D  T5=260550.  *(FSS(IC+5)+FSS(  K+7)) 

D  T6*427368.  *(FSS(IC+6)  ) 

D  ATI* (5**0»00l ) V299376* 

D  AT2-T1+T2-T3+T4-T5+T6 

D  A*ATl*AT2 

L*1002-I 
M*1002-I+l 
D  250  GS(L  V*GS(M)+A 


J 
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DO  260  1*11*1001 
J*  1001-1 

T  1*16067 •  *(FS  (J+D+FS  U+ll)) 

T2*1Q6300#  * ( FS  ( J+2)+FS  (J+10)) 

T3*48525.  *(FS  (J+31+FS  (J+9>) 

T4*272400.  #(FS  (J+4J+FS  ( J+8 > > 

T5*260550.  * ( FS  (J+5J+FS  (J+7M 
T6»427368.  *(FS  (J+6)) 

ATI* ( 5»*0*01  )/299376. 

AT2*T1+T2-T3+T4~T5+T6 
A*AT1#AT2 
K=1002-I 
L=1002-I+10 
260  GS(K)*GS(L)+A 

-  COMPUTE  FLS 

DO  270  1*1*1001 
270  FLS<n*FSm-TC  I  )*GSU) 

-  PRINT  DISTRIBUTIONS 

WRITE  OUTPUT  TAPE  6*51*NU 

WRITE  OUTPUT  TAPE  6*15* (FS  (11*1*1*1001*100) 

WRITE  OUTPUT  TAPE  6*52»NU 

WRITE  OUTPUT  TAPE  6.15»(GS  ( I ) » 1*1 .1001 *100 ) 

WRITE  OUTPUT  TAPE  6*53. NU 

WRITE  OUTPUT  TAPE  6*15* (FLS( I) *1*1*1001*100) 

-  READ' MAXIT 

READ  INPUT  TAPE  5*10»MAXIT 
WRITE  OUTPUT  TAPE  6*10»MAXIT 

-  READ  INDEX  IOFD  INDEX  ON  D 

READ  INPUT  TAPE  5. 10* IOFD1 . IOFD2 » IOFD3 
WRITE  OUTPUT  TAPE  6*10* IOFD1 ♦ IOFD2. IOFD3 

DO  2000  IOFD  *  IOFD1 ♦ IOFD2 ♦ IOFD3 

-  READ  AND  PRINT  THCT ( IOFD) *  I0FLl*I0FL2* IOFL3 

READ  INPUT  TAPE  5*15*THCT( IOFD) 

READ  INPUT  TAPE  5. 10* I0FL1 ,» I0FL2 » I0FL3 

WRITE  OUTPUT  TAPE  6*21 

WRITE  OUTPUT  TAPE  6.15 *THCT ( IOFD) 

WRITE  OUTPUT  TAPE  6.10* I0FL1 » I0FL2* I0FL3 

-  COMPUTE  ALAM  *  F(ATH)  WHERE  ATH  *  OR  GREATER  THAN  THCT (IOFD) 

-  NOTE  ATH  IS  STORED  IN  DECREASING  ORDER 

DO  400  J*1 *360 

I F ( ATH( J )-  THCT( IOFD) )350*350»370 
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350  JLONG*=J-l 
GO  TO  425 
C 

370  DTH=D( IOFD) /ATH ( J ) 

HD*  XLOCF  (L(  ID 
L ( 2 ) *  XLOCF (  T  ( 1 ) ) 

L ( 3 ) *  XLOCF I FS! 1 ) ) 

L  ( 4 )  *  1 
H  5  )  =  1 
L  ( 6 )  =  5 
L  ( 7  )  =  1001 
L  ( 8  )  =  0 

FSDTH*TABF(DTH*L( 1 ) ) 

I  — L  C  8 ) 

GO  TO! 390*375 ) » I 
C 

375  PRINT  380*  IOFD*  J*  DTH 

380  FORMAT! 15H1TFS! D/ATH)  BAD2 1 5  *E20*8 ) 

GO  TO  9999 
C 

390  TERM1*1. 

TERM2* ( 2 •*  ATH ( J ) )  /  (((  l.-ATHI J)**2  )**2  )*  FSDTH ) 
ALAM(J)*  TERM1  *  TERM2 
400  CONTINUE 
JLONG=360 
C 

C  -  COMPUTE  CTH  *  F ( CLAM)  *  WHERE 

C  -  CLAM  *  10**L  »  BY  USING  ALAM*  F(ATH)  TABLE 

C 

425  DO  750  U=I0FL1*I0FL2*I0FL3 
CLIJ-I*  BLIJ) 

CLAM ( J ) *  BLAM(J) 

L  ( 1  >  *  XLOCF!  LID) 

L(2>  =  XLOCF ( ALAMI 1 ) ) 

L ( 3 ) *  XLOCF ( ATH! 1 ) ) 

L  <  4 )  *1 
L ! 5 ) =1 
H  6 )  “5 
LI?)*  JLONG 
L  (  8 ) *0 

CTH I J ) =  T ABF t CLAM (J) *L(  ID 
GO  TO  (750*710  )*I 
710  PRINT  720 »  IOFD*  J»CTH(J) 

720  FORMAT! 16H1CTH»F(CLAM)  BAD2 1 5 *E20.8 ) 

GO  TO  9999  . 

C 

750  CONTINUE 
C 

C  COMPUTE  CRHOT 

C 

DO  800  J»10FL1»IOFL2»IOFL3 

CRHOT ( J ) «  CTH! J)**2  /  (1.-  CTH!J)**2  ) 

800  CONTINUE 
C 

C  COMPUTE  CRHO  BY  NEWTON  RAPHSON 

C 

DO  900  J*IOFLl*IOFL2»IOFL3 
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RHONEW-CRHOT < J) 

AGPRN-O. 

C 

DO  860  IT-1 * MAX IT 
IF< IT-10)825, 825*830 
825  RHO-RHONEW 
GO  TO  850 
830  RHO-RHONEW 

I F ( AGPRO-AGPRN ) 835 t 835 »84Q 
835  GO  TO  890 
840  GO  TO  850 
C 

850  AGPRO-AGPRN 
V*D( IOFD) 

FLAM-10. **CLU) 

THETA«SQRTF(RH0/(1.+RH0) ) 

CON-1. /SQRTF<2. *3. 141592653589793)  ' 
DTH-V/THETA 

FSDTH*CON*EXPF( (-.5 )*DTH**2 ) 

T1=C.5*FLAM)  /  (SQRTF (RHO»CRHO+l.)**3) ) 
T2-1 . 

GPR-(T1*T2*FSDTH)-1. 

TA-FLAM/4. 

TB1»S0RTF(RH0*(RH0+1.) ) 

TB2-TB1**5 

TB-1./TB2 

TQ-( (-4. )*RH0**2 )  +  (V**2-l.)*RHO  +  V**2 
G2PR-T  A*TB*FSDTH*TQ 
RHONEW-RHO  -  (GPR/G2PR ) 

AGPRN-ABSF (GPR ) 

C 

880  CONTINUE 
C 

890  CRHO ( J) -RHONEW 
ITER( J) -IT 
CGPR ( J ) -GPR 
C 

900  CONTINUE 

-  COMPUTE  CLGRHO 

DO  950  J*IOFLl»IOFL2*lQFL3 
IF (CRHO(J)-O. >910*910*920 
910  CLGRHO!  IOFD* J)  — 99.9999 
GO  TO  950  , 

920  CLGRHO ( IOFD* J ) -LOGIOF (CRHO( J )) 

950  CONTINUE 

-  COMPUTE  CNG*  CLGNG 

1000  DO  1100  J-  IOFL1 ♦ 10FL2 ♦ IOFL3 
DTH»D(I0FD)/CTH(J> 

HD*  XLOCF (  Kill 
L(2)»  XLOCF  (Til)) 

L ( 3 ) *  XLOCF ( FLS ( 1 ) ) 

L 1 4 )  *  1 
L  ( 5  » -  1 
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L  ( 6 )  *  5 
L(7>«  1001 
L(8)»  0 

flsdth=tabf(dth»u  m 

I *L( 8 ) 

60  TO  (1050, 1020), I 
1020  PRINT  1030,  I0FD,  J,  FLSDTH 
1030  FORMAT ( 15H1FLSDTH  BAD2 1 5 ,E20.8 ) 

60  TO  9999 
C 

1050  CN6(J>*  CLAM(J)  *  CTH(J)  *  FLSDTH  -  CRHO(J) 

IF(CNG(J>-0, ) 1060,1060*1070 
1060  CL6N6 ( 1 0FD , J 1  *  -99.9999 
60  TO  1100 

1070  CL6NG ( IOFD, J ) »  L0610F  (CNG(J)J 
C 

1100  CONTINUE 

-  RHO  ASYMPTOTIC  CALCULATION  FOR  VARIOUS  CLAM  VALUES 

KOFD*( 10*I0FD)-9 
FSOFD-FS(KOFD) 

DO  1600  J=IOFL1,IOFL2*IOFL3 
1600  CRHOAS( J)=SQRTF( .5*CLAM( J)*FS0FD) 

DO  1650  J=IOFLl*IOFL2»IOFL3 
1650  RHORAT ( J ) =CRHO ( J ) /CRHOAS ( J ) 

-  WRITE  FUNCTIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6»81*D(I0FD) 

WRITE  OUTPUT  TAPE  6,61 
WRITE  OUTPUT  TAPE  6,22 

WRITE  OUTPUT  TAPE  6*71, ( CL (J) »CLAM( J) *CTH( J ) «CRHO( J ) »CLGRHO  < IOFD ♦ 
XJ ) ,CNG( J ) ,CLGNG ( IOFD* J ) * J* I0FL1 , IOFL2, IOFL3 ) 

-  WRITE  CLGRHO  AND  CLGNG  ON  TAPE  2 

WRITE  OUTPUT  TAPE  2 ,15, ( CLGRHO< IOFD, J) , J»I0FL1, IOFL2 , I0FL3 ) 

WRITE  OUTPUT  TAPE  2*15*(CLGNG  ( IOFD, J) , J* IOFL1 ♦ IOFL2 , I0FL3 ) 

-  WRITE  CRHO  FUNCTIONS 

WRITE  OUTPUT  TAPE  6,81*D(I0FD) 

WRITE  OUTPUT  TAPE  6,22 

WRITE  OUTPUT  TAPE  6,72  * < CL( J ) ,CRH0T ( J) ,CRHO( J ) ,CGPR ( J ) , I  TER ( J ) , J* I 
XOFL1 ♦ IOFL2 » IOFL3 ) 

-  WRITE  ASYMPTOTIC  CALCULATIONS  ON  TAPE  6 

WRITE  OUTPUT  TAPE  6,81*D(I0FD) 

WRITE  OUTPUT  TAPE  6*63 
WRITE  OUTPUT  TAPE  6,22 

WRITE  OUTPUT  TAPE  6 ,73 ,  ( CL ( J)  *CRHO(  J)  ,CRHOAS<  J )  , RHORAT ( J ) « J*>  IOFL 1 
X , IOFL2, IOFL3 ) 

- COMPUTE  PHI  AND  PS  I  FOR  D 
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C 

PHIAS( IOFD)«.5*LOG10F(.5*FSOFD) 

C 

FLSOFD*FLS(KOFD> 

PSIAS( IOFD) *L0G10F (FLSOFO) 

-  INDEX  ON  D 

2000  CONTINUE 

-  WRITE  SS  AND  NG*  PLUS  PHIAS  AND  PSIAS 

CALL  WRSSNG ( NU  »  CL » CLGRHO  *  CLGNG ♦ PH  I  AS  * PS  I AS ) 

- END 

9999  CONTINUE 

END  FILE  2 
END  FILE  8 
REWIND  2 
REWIND  8 
PRINT  6 
PAUSE 
CALL  EXIT 
END 
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D.4  Subprogram  for  Printing  R  and  T. 

-  OPTIMAL  SAMPLE  SIZE  AND  NET  GAIN  EDIT  AND  PRINT 

SUBROUTINE  WRSSNG  <NU»CL *CLGRHO»CLGNG*PHIAS»PSlAS) 

DIMENSION  CL ( 121 ) yCLGRHO (41*121) *CLGNG( 41*121) 

COMMON  CL*CLGRHO*CLGNG 
DIMENSION  PHIAS! 121 ) *PSI AS! 121 ) 

COMMON  PHIAS*PSIAS 

16  FORMAT! 3110 ) 

36  FORMAT ! 1H  ) 

84  FORMAT !F10*1*2X*11F9»4) 

85  FORMAT ! 1H16X *2HSSI 3  *34X» 19H0PTIMAL  SAMPLE  SIZE) 

86  FORMAT! 1H150X.8HNET  GAIN48X»2HNGI3) 

87  FORMAT! 1H016X*3H0.06X*3H0.26X*3H0.46X*3H0*66X»3H0.86X»3H1.06X»3H1. 
X26X*3H1.46X*3H1.66X»3H1*86X*3H2*0) 

88  FORMAT! 1H016X»3H2.06X *3H2*26X»3H2.46X»3H2 .66X .3H2.86X.3H3.06X *3H3 . 
X26X*3H3.46X*3H3.66X*3H3.86X*3H4.0) 

95  FORMAT! 10H  PHI2X.1 IF9.4) 

96  FORMAT! 10H  PSI2X.11F9.4) 

-  WRITE  SS  AND  NG 

READ  INPUT  TAPE  5 ♦ 16.LP1 »LP2 *LP3 
WRITE  OUTPUT  TAPE  6.16*LP1*LP2.LP3 

-  SS  FIRST  PAGE 

WRITE  OUTPUT  TAPE  8*  85 ♦  NU 

WRITE  OUTPUT  TAPE  8.  87 

WRITE' OUTPUT  TAPE  8»  36 

C 

DO  2200  I«LP1»LP2»LP3 
C 

WRITE  OUTPUT  TAPE  8 *84* CL ( I ) » (CLGRHOI J *  I ) *  J«l»21*2  ) 

C 


IF! 

1-21)2200* 

2150* 

2026 

2026 

IF! 

1-26)2200* 

2150* 

2031 

2031 

IF! 

1-31)2200* 

2150* 

2036 

2036 

IF! 

1-36)2200* 

2150* 

2041 

2041 

IF! 

1-41)2200* 

2150* 

2046 

2046 

IF! 

1-46)2200* 

2150* 

2051 

2051 

IF! 

1-51)2200* 

2150. 

2056 

2061 

IF! 

1-61)2200* 

2150* 

2066 

2066 

IF! 

1-66)2200* 

2150* 

2071 

2071 

IF! 

1-71)2200* 

2150* 

2076 

2076 

IF! 

1-76)2200* 

2150* 

2081 

2081 

IF! 

1-81)2200* 

2150* 

2086 

2086 

IF! 

1-86)2200* 

2150* 

2091 

2091 

IF! 

1-91)2200* 

2150* 

2096 

2096 

IF! 

1-96)2200* 

2150* 

2200 

2150  WRITE  OUTPUT  TAPE  8*  36 

C 

2200  CONTINUE 
C 

WRITE  OUTPUT  TAPE  8*  36 
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WRITE  OUTPUT  TAPE  8*  36 

WRITE  OUTPUT  TAPE  8 ,95 , < PHI  AS < J ) , J= 1 ,21  ,2 ) 

-  NG  FIRST  PAGE 

WRITE  OUTPUT  TAPE  8,  86,  NU 

WRITE  OUTPUT  TAPE  8,  87 

WRITE  OUTPUT  TAPE  8,  36 

DO  2400  I =LP 1 »LP2  »  LP3 

WRITE  OUTPUT  TAPE  8,  84,  CLU)  •  (CLGNGI J  » I )  *  J*1  ,21 ,2  > 


I F  ( 

1-21)2400, 

2350* 

2226 

2226 

IF  ( 

1-26)2400, 

2350  * 

2231 

2231 

I F  ( 

1-31 J2400, 

2350, 

2236 

2236 

I F  ( 

1-36)2400, 

2350, 

2241 

2241 

I  F  ( 

1-41)2400, 

2350, 

2246 

2246 

I F  ( 

1-46)2400, 

2350, 

2251 

2251 

I  F  ( 

1-51)2400, 

2350* 

2256 

2256 

I  F  ( 

1-56)2400, 

2350, 

2261 

2261 

I F  ( 

1-61 ) 2400  » 

2350, 

2266 

2266 

I F  ( 

1-66)2400, 

2350, 

2271 

2271 

I  F  ( 

1-71)2400, 

2350, 

2276 

2276 

I F  ( 

1-76)2400, 

2350, 

2281 

2281 

I F  ( 

1-81)2400, 

2350, 

2286 

2286 

I F  C 

1-86)2400, 

2350, 

2291 

2291 

I F  ( 

1-91)2400, 

2350, 

2296 

2296 

I F  ( 

1-96)2400, 

2350, 

2400 

2350  WRITE  OUTPUT  TAPE  8,  36 


2400  CONTINUE 

WRITE  OUTPUT  TAPE  8,  36 

WRITE  OUTPUT  TAPE  8,  36 

WRITE  OUTPUT  TAPE  8 ,96, ( PSI AS< J ) » J*1 ,21 ,2 ) 

-  WRITE  SS  AND  NG»  SECOND  PAGE 

READ  INPUT  TAPE5»16»LPl,LP2»LP3 
WRITE  OUTPUT  TAPE  6 , 16.LP1  ,LP2 ,LP3 

-  SS  SECOND  PAGE 

WRITE  OUTPUT  TAPE  8,  85,  NU 

WRITE  OUTPUT  TAPE  8,88 

WRITE  OUTPUT  TAPE  8,  36 

C 

DO  2600  I»  LP1 *LP2,LP3 
C 

WRITE  OUTPUT  TAPE  8,  84»CL( I ) » (CLGRHOI J, I ) ,  J*21.41,2> 

1  C 

IF  (I-  21)2600,  2550,  2426 
2426  IF  (I-  26)2600,  2550,  2431 
2431  IF  (I-  31)2600,  2550,  2436 
2436  IF  (I-  36)2600,  2550,  2441 
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2441 

IF 

( I-  41)2600, 

2550. 

2446 

2446 

IF 

l I-  46)2600. 

2550. 

2451 

2451 

IF 

( I-  51)2600. 

2550. 

2456 

2456 

IF 

(I-  56)2600, 

2550. 

2461 

2461 

IF 

(I-  61)2600, 

2550. 

2466 

2466 

IF 

( I-  66)2600, 

2550. 

2471 

2471 

IF 

(I-  71)2600, 

2550. 

2476 

2476 

IF 

( I-  76)2600. 

2550, 

2481 

2481 

IF 

( I-  81)2600, 

2550. 

2486 

2486 

IF 

( I-  86)2600, 

2550. 

2491 

2491 

IF 

( I-  91)2600, 

2550. 

2496 

2496 

IF 

(I-  96)2600, 

2550. 

2501 

2501 

IF  ( 

1-101 )2600, 2550, 2506 

2506  I F ( I  — 106 ) 2600 *2550*2511 
2511  IF( I  — 1 1 1 )2 600,2550,2516 
2516  I F (1-116 ) 2600 » 2550  *2600 

< 

2550  WRITE  OUTPUT  TAPE  8»  36 
C 

2600  CONTINUE 
C 


WRITE  OUTPUT  TAPE  8*  36 

WRITE  OUTPUT  TAPE  8,  36 

WRITE  OUTPUT  TAPE  8»95»<PHIAS(J1 »J321»41»2) 


-  NG  SECOND  PAGE 

WRITE  OUTPUT  TAPE  8.  86.  NU 
WRITE  OUTPUT  TAPE  8»88 

WRITE  OUTPUT  TAPE  8 ,  36 

C 

DO  2800  1=  LPl.LP2.LP3 

C 

WRITE  OUTPUT  TAPE  8*  84.  CL ( I ) . ( CLGNG< J . I ) . J»21 ,41 *2 ) 
C 


IF 

( 

I- 

21) 

2800  * 

2750. 

2626 

2626 

IF 

( 

I- 

26) 

2800  , 

2750, 

2631 

2631 

IF 

( 

I- 

31) 

2800  , 

2750. 

2636 

2636 

IF 

( 

I- 

36) 

2800  , 

2750, 

2641 

2641 

IF 

( 

I- 

41)2800, 

2750, 

2646 

2646 

IF 

( 

I  - 

46) 

2800  • 

2750, 

2651 

2651 

IF 

( 

I  ~ 

51) 

2800  , 

2750, 

2656 

2656 

IF 

( 

I- 

56) 

2800, 

2750, 

2661 

2661 

IF 

( 

I- 

61) 

2800, 

2750, 

2666 

2666 

IF 

( 

I- 

66) 

2800, 

2750, 

2671 

2671 

IF 

( 

I- 

71) 

2800, 

2750, 

2676 

2676 

IF 

( 

I- 

76) 

2800, 

2750, 

2681 

2681 

IF 

( 

I- 

81) 

2800, 

2750, 

2686 

2686 

IF 

( 

I- 

86) 

2800, 

2750, 

2691 

2691 

IF 

{ 

I~ 

91) 

2800, 

2750, 

2696 

2696 

IF 

( 

I- 

96) 

2800, 

2750, 

2701 

2701  I F ( 1-101)2800.2750*2706 
2706  I F ( 1-106)2800.2750.2711 
2711  I F ( 1-111)2800.2750.2716 
2716  IF( 1-116) 2800. 2750.2800 
C 

2750  WRITE  OUTPUT  TAPE  8,  36 
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C 

2800  CONTINUE 
C 

WRITE  OUTPUT  TAPE  8.  36 
WRITE  OUTPUT  TAPE  8*  36 

WRITE  OUTPUT  TAPE  8 .96* ( PSI AS< J ) » J=21 »A1 »2 ) 

- END 

RETURN 
END 
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C 

C 

10 

25 

30 

C 

C 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

C 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

D 

C 

D 

D 

D 


A  and  R  for  v  <  « 
c _ c _ 


CUTOFFS*  NU  LESS  THAN  INF 


DIMENSION 

DIMENSION 

DIMENSION 


T( 1001),  FS(lOOl) »GS( 1001) *FLS( 1001) 
TS(91) *FSS(91) 

FLNUT ( 19 ) ,STVT ( 19 ) ,BETAT ( 19 ) 


DIMENSION  L( 11 ) 

DIMENSION  D( 101 ) ,THCTT ( 101 ) * THCT ( 101 ) tFLAMCT ( 101 ) »FLCT ( 101 ) 
DIMENSION  RHOCT(lOl) *RCT(101) 


DIMENSION  TH1  ( 100.1 )  * 


Dl(lOOl) 


FORMAT (7110) 
FORMAT ( 5E20*8 ) 
FORMAT ( 1H1 ) 

SPECIFY  CONSTANTS 


FLNUT ( 3 )=2* 
FLNUT ( 4 ) *3 • 
FLNUT ( 5 ) *4* 
FLNUT ( 6  )*5* 
FLNUT  ( T )  *=6  • 
FLNUT ( 8  > “7 • 
FLNUT (9 )=8* 
FLNUT ( 10 ) =9* 
FLNUT(11)=10. 
FLNUT < 12 )  =  12* 
FLNUT (13 )=15* 
FLNUT(14)=20. 
FLNUT"!  15. )  *24* 
FLNUT (16 )=30. 
FLNUT ( 17 )=40* 
FLNUT(18)=60. 
FLNUT ( 19 )=120* 


STVT 1 3 )  = 
STVT ( 4 ) = 
STVT 1 5 )  = 
STVT (6) = 
STVT ( 7 )  « 
STVT (8)* 
STVT (9 ) = 
STVT (10) 
STVT (11) 
STVT ( 12 ) 
STVT (13) 
STVT (14) 
STVT (15) 
STVT (16) 
STVT ( 17 ) 
STVT (18) 
STVT (19) 


0* 492622851 1662845 E-02 
Q« 1064199529207075 E-02 
0*28100181 13579956E-03 
0* 8547 378787 148 180E-04 
0* 2895991 3774768 13 E-04 
0* 1069710 1445 38641E-04 
0* 4244090 76 38 14246 E-05 
=*0*  1789 11 87 15962368  E-05 
=0 *7947 76 5877982 060 E— 06 
*0. 1790661843836187E-06 
*0 *249844907 14624 72 E—07 
=0* 158 18908 79 3 57 194 E~08 
=0* 2457 7627638 3765 3 E“09 
=0.22876257041 14807E-10 
=0.9656558502057788E-12 
=0* 1068 84295 73 3 7 04 IE— 13 
=0.85697205876S0508E-17 


BETAT ( 3 ) =  0* 200000000000000 0E+01 
BETAT ( 4 ) =  0.1570796326794897E+01 
BETAT ( 5 ) *  0*1333333 33 3333333E+01 
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BETAT ( 6 )  =  0,1178097245 09 6172E+01 
BETAT { 7 )  =  0, 1066666666666667E+01 
BETAT ( 8 ) =  0 ,98 174 77 042 4681 04 E+00 
BETAT ( 9 )  =  0,9142857 1428 5714 3E+00 
BETAT ( 10) =0,8590292412 159591E+00 
BETAT ( 11 ) =0.8 126984 126984 12 7 E+00 
BETAT ( 12) =0.7388167388167388E+00 
BETAT ( 13) =0.6580777580029401E+00 
BETAT ( 14) =0,567 546385 5 03 041 9 E+00 
BETAT (15 )=0, 5170194816 176779E+00 
BETAT ( 16 ) =0 *4614 7455 34009741 E+00 
BETAT ( 17 ) =0.3988 1730689488 10E+00 
BETAT ( 18  >=0,32 49 5 542 03948302 E+00 
BETAT ( 19 )=0,2293000137934539E+00 

-  GENERATE  T  AND  TS 


I  T  =  0 
T  ( 1 )  =  I T 

DO  110  1=2*1001 
I T=  IT  +  1 
TT=IT 

110  T ( I )  =  TT/100. 

IT  =  9910 
TS( 1 )=  IT/1000 
DO  115  1=2*91 

I T= I T+l 
TT=IT 

115  TS< I)=TT/1000. 

- D  LIST  OUlllO 

ID  =  0 
D ( 1 ) =  ID 
DO  130  1=2,101 
ID=ID+1 
TD=  ID 

130  D ( I ) =  TD/10 • 

-  READ  NUMNUS  SET  UP  MAIN  INDEX  IMAIN 

READ  INPUT  TAPE  5»10»NUMNUS 
WRITE  OUTPUT  TAPE  6»10*NUMNUS 

DO  4000  IMAIN=1 »NUMNUS 


- READ  INU »NU, I0FD1 » I0FD2 , I0FD3 

READ  INPUT  TAPE  5* 10, INU»NU»I0FD1 ,I0FD2 ,I0FD3 
WRITE  OUTPUT  TAPE  6*30 

WRITE  OUTPUT  TAPE  6 ,10* 1NU*NU* I0FD1 * I0FD2 . I0FD3 

FLNU*FLNUT( INU) 

STV=STVT( INU) 

BETA=BETAT ( INU ) 
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C 

C  -  COMPUTE  FS 

C 

D  TERM1*  1  •  /  ( BETA*SQRTF (FLNU) ) 

D  EXP*~( ( FLNU+1 • ) /2 • ) 

C 

DO  230  1*1.1001, 

D  230  FS  ( I )  =TERM1*  ((l.+(T  (I )**2 > /FLNU) **EXP ) 

C 

C  -  COMPUTE  GS 

C 

D  GS ( 1001 )  *  STV 

C 

DO  235  1*1.91 

D  235  FSS( I )=TERM1*  ( ( 1 .+< TS< I )**2 > /FLNU)**EXP ) 

C 

DO  250  1=2.10 
J*  11-1 
K*  10*<J-1) 

D  Tl*16067.  *(FSS(K+1)+FSSU+11> ) 

D  T2*106300.  *{ FSS( K+2 )+FSS( K+10 ) ) 

D  T3=48525.  *(F$S(K+3)+FSS(K+9) ) 

D  T 4=272400*  *(FSSC K+4)+FSS(K+8)) 

D  T5*260550.  * l FSSOC+5 )+FSS( K+7 ) ) 

D  T6=427368.  *(FSS(K+6)) 

D  ATI* ( 5»*0»001 ) Z299376 • 

D  AT2=T1+T2-T3+T4-T5+T6 

D  A*AT1*AT2 

L*1002-I 
M*1002-l+l 
D  250  GS(L)*GS(M)+A 
C 

DO  260  1*11.1001 
J*  1001-1 

D  Tl  =  16067.  *<FS  (J+ll+FS  (J+ll>> 

D  T2*106300.  * ( FS  (J+2)+FS  (J+10M 

D  T3*48525.  *(FS  <J+3)+FS  (J+9)) 

D  T4*272400.  *<FS  <J+4)+FS  ( J+8 ) ) 

D  T5*260550.  *<FS  (J+5)+FS  (J+7)> 

D  T6*427368.  *CFS  (J+6)) 

D  AT1=<5.*0.01  )/299376. 

D  AT2*T1+T2-T3+T  4-T5+T6 

D  A*AT1#AT2 

K*1002-I 
L*1002-I+10 
D  260  GS(K)*GS(L)+A 
C 

C  —  COMPUTE  FLS 
C 

D  BOT  *  FlNU-1 • 

DO  270  1*1.1001 

D  270  FLS( I ) = ( (FLNU+  T ( I )**2. ) /BOT )*FS<I )  -  T( 1 )*  GS(I) 
C 

C  -  WRITE  FS.GS.FLS 

C 

WRITE  OUTPUT  TAPE  6.25. ( FS( 1> » 1*1 .1001*100) 

WRITE  OUTPUT  TAPE  6.25 .( GSU  ).  1*1  *1001 . 100) 
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WRITE  OUTPUT  TAPE  6 .25 » (FLS ( I ) *  I -1 .1001  *  100 ) 

C 

C - ROUTINE  FOR  COMPUTING  THCTT  *THCT  .FLAMCT  »FLCT  *RHOCT  »AND  RCT 

C 

C  -  DECISION  FOR  NU  ROUTINE 

C 

I F ( NU— 2 )290 » 290*300 
C 

C - NU  *  2 

C 

290  DO  295  I =I0FD1 * IOFD2 * IOFD3 
THCT ( I  )  =  0* 

FLAMCT ( I )*  2.*  D(  I  ) 

IF  (FLAMCTU)-O.  >293*293*291 

291  FLCT(I)  =  LOGIOF ( FLAMCT ( I )) 

GO  TO  295 

293  FLCT ( I )  =  -• 99999999 E+09 
295  RCT ( I ) =-0*99999999E+09 
GO  TO  3000 

-  COMPUTE  TH1  AND  ASSOCIATED  Dl 

300  ASSIGN  310  TO  ISW 
DO  330  1=1*1001 

PHI  *  T(I)  *  ( (FLNU-1.)/(FLNU  +  T<I>**2))  *  <GS( I >/FS< I ) ) 

I F (  PHI-  .5)  330*330*305 
305  GO  TO  ISW.  (310.315) 

310  ASSIGN  315  TO  ISW 
I  FIRST  =  I 
LENGTH* 1002-1  FIRST 
315  THl(l)  *  SORTF  <2.*PHI-  1.) 

Dl( I)  =  TH1 ( I )  *  T<I> 

330  CONTINUE 

WRITE  OUTPUT  TAPE  6*10. IFIRST .LENGTH 
WRITE  OUTPUT  TAPE  6.  25  *  *  Dl ( I ) » 1*1 *1001 ) 

WRITE  OUTPUT  TAPE  6*  25, ( THl ( I) . 1*1  *  1001 ) 

- COMPUTE  THCTT  BY  INTERPOLATION  IN  Dl*  F(TH1)  TABLE  * 

DO  450  J*IOFDl,IOFD2»IOFD3 
IF(J— 1)350*350*360 
350  THCTT ( J ) =0* 

GO  TO  425 
360  DT*D(J) 

DO  390  1=1*1001 
IF(DT-D1(I )) 370.380.390 
370  TA»DT-D1(I-1) 

TB=D1 ( I )-Dl( I-l) 

THCTT(J)-THKI-l)  +  (TA/TB)*(TH1(I)-TH1(I-U  ) 

GO  TO  425 

380  THCTT ( J ) *THl (  I ) 

GO  TO  425 
390  CONTINUE 

THCTT ( J ) *THl < 1001 ) 

GO  TO  425 
C 
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425  WRITE  OUTPUT  TAPE  6*25 » THCTT ( J ) 

C 

450  CONTINUE 

-  USING  THCTT  COMPUTE  THCT  BY  NEWTON-RAPHSON 


C 

1100 

1110 


C 

1120 

1125 


C 

1140 


DO  1500  J=  I0FD1*I0FD2»10FD3 
I  ND*0 

T I =D ( J ) /THCTT ( J) 

WRITE  OUTPUT  TAPE  6.11l0*TI 
FORMAT ( E20.8 ) 

L  ( 1 )  *  XLOCF(Ul)  ) 

L  ( 2  )  a  XLOCF ( T ( 1 ) ) 

L  (  3 )  =  XLOCF(GSU)) 

L  ( 4  >  =  1 
L  (  5 )  =  1 
L  ( 6 )  =  5 
L ( 7 ) =  1001 
L  (  8  )  *  0 

GOFTI  *  TABF (  TI  *L(1) > 

I aL( 8 ) 

GO  TO  I  1140*1120) »1 

PRINT  1125*  TI 
FORMAT (  6H1G  BADE20.8) 

GO  TO  9999 

L(l)-  XLOCF (L ( 1 )  ) 

L  1 2 ) -  XLOCF ( T  ( 1 ) ) 

L ( 3 ) -  XLOCF ( FS ( 1 )  ) 

L  (  4  )  =  1 

L  ( 5 )  a  1 
L  ( 6 )  =  5 

L(7J-  1001 

L  (  8 )  a  0 

FOFT I =  TABF ( TI *L ( 1 )  ) 

I«  L ( 8 ) 

GO  TO  (1180*1160) »I 


C 

1160  PRINT  1165*  TI 

1165  FORMAT  ( 6H1F  BADE20.8) 

GO  TO  9999 
C 

1180  0 1  =  T I  * ( ( FLNU-1 • ) / ( FLNU+  TI**2))  *  (GOFTI  /FOFT I ) 

Yl  =  QI  -  • 5* ( 1 •  +  (D(J)/TI )**2) 

QIPR»( (FLNU-1. )/(FLNU+  TI#*2) )  *< ( FLNU/ ( FLNU+T 1**2 ) ) * ( l.+T I **2 ) * 
X  ( GOFTI /FOFT I )  -  TI  ) 

YIPR=  QIPR+  (  D(J)«*2  /  T 1**3 ) 

T INC®  TI  -  (YI/YIPR) 


C 

IF ( IND-O) 1200* 1200* 1220 
1200  IND*  1 

TI  ■  TINC 
GO  TO  1100 
C 

1220  I F (  TINC-  TI ) 1300*1300*1250 
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1250 

T I =  TINC 

GO  TO  1100 

1300 

THCTU)  =  D(J )  /T I 
WRITE  OUTPUT  TAPE 

6*1305 

1305 

FORMAT (1H  ) 

WRITE  OUTPUT  TAPE 

6  *1310  »THCT ( J ) 

1310 

FORMAT ( E20 • 8 ) 

1500 

CONTINUE 

—  COMPUTE  FLAMCT*  FLCT  USING  FS  TABLE 

DO  2500  J=  I0FD1 • 10FD2* IOFD 
DTH  *  D(J)  i  THCTU) 

L ( 1 ) *XLOCF ( L ( 1 ) ) 

L ( 2 ) “XLOCF ( T  ( 1 ) ) 

L  (  3 ) =XL0CF ( FS( 1 ) ) 

L<4)=1 
L  ( 5  )  =1 
L  ( 6 )  =5 
L  ( 7 )  *  1001 
L  ( 8 )  *  0 

FSDTH  *  TABF(DTH»L<1) ) 

1=  L  ( 8 ) 

GO  TO  (2440*2420) *  I 
2420  PRINT2425*J*DTH 

2425  FORMAT  (  15H1FS( D/THCT )  BADI4.E20.8) 

GO  TO  9999 

2440  FLAMCA* ( FLNU-1 • ) / ( FLNU+DTH**2 ) 

FLAMCB  *  <2.*THCT(J))/  ( ( ( l.-THCT( J)**2  )**2  )*  FSDTH) 
FLAMCT ( J ) =FLAMCA*FLAMCB 
IF (FLAMCT(J)-O.) 2450*2450*2445 
2445  FLCT ( J)  =  LOG10F  ( FLAMCT ( J ) ) 

GO  TO  2500 

2450  FLCT(J)  *  -.99999999E+09 
2500  CONTINUE 

-  COMPUTE  RHO  AND  L06RH0 

DO  2950  J=I0FD1 » I0FD2  *IOFD3 

RHOCT U )*THCT ( J)**2  / ( 1 .-THCT ( J) **2  ) 

IF (RHOCT(J)-O. >2920*2920*2930 
2920  RCT(J) *-99.9999 
GO  TO  2950 

2930  RCT ( J ) *LOG10F ( RHOCT < J ) ) 

2950  CONTINUE 

- PRINT  D *  THCTT*THCT*FLCT.RCT 

3000  WRITE  OUTPUT  TAPE  6»3050*NU 
3050  FORMAT ( 4H1NU* 1 4 ) 

WRITE. OUTPUT  TAPE  6*3100 


-  62  - 


non  non  non 


Program  E. 1 ,  page 


3100  FORMAT ( 110H0 

D 

THCTT 

X  THCT 

FLCT 

’  ^  l 

RCT ) 

DO  3200  J«I0FD1*I0FD2*I0F03 

WRITE  OUTPUT  TAPE  6*3150 *D< J ) *THCTT ( J) *THCT ( J ) .FLCT ( J ) *RCT ( J ) 
3150  FORMAT<F10.1»4E25.8) 

3200  CONTINUE 

—  Blank  out  di  and  thi 

DO  3500  J=l*1001 
Dl ( J  ) *0 • 

TH1(J)*0. 

3500  CONTINUE 

-  MAIN  INDEX 

4000  CONTINUE 

- END 

9999  CONTINUE 
CALL  EXIT 
END 
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E.2  A  and  R  for  v  =  » 
c _ c _ 


-  CUTOFFS*  NU  = I NF 

DIMENSION  T ( 1001 ) *  FS ( 1001 ) *GS ( 1001 ) *FLS ( 100 1  ) 

DIMENSION  T S ( 9 1 ) »FSS( 91 ) 

DIMENSION  L( 11) 

DIMENSION  D (  101 )  •  THCTT  ( 101 )  »THCT(10l)  »FLAMCT<101 >  *FLCT(101 ) 
DIMENSION  RHOCT ( 101 ) *RCT  < 10 1 ) 

DIMENSION  THl(lOOl).  DK1001) 

10  FORMAT (7110) 

25  FORMATt 5E20.8  ) 

30  FORMAT ( 1H1 ) 

-  NU-999 

NU=999 

-  READ  AND  PRINT  IOFD 

READ  INPUT  TAPE  5 » 10* I0FD1 » I0FD2* IOFD3 
WRITE  OUTPUT  TAPE  6*10* IOFDl*IOFD2* IOFD3 

-  GENERATE  T  AND  TS 


I  T=0 
T  ( 1 )  =  I T 

DO  110  I =2  *  1001 
IT=  IT  +  1 
TT=IT- 

110  Ttl)  *  TT/100. 

IT  *  9910 
TS(1)«  IT/1000 
DO  115  1=2*91 

I T= I T+l 
TT=IT 

115  TS( I )=TT/1000. 

-  D  LIST  01.1)10 

ID  =  0 
D ( 1 ) =  ID 
DO  130  1=2,101 
ID=ID+1 
TD=  ID 

130  D ( I ) *  TD/10. 

-  COMPUTE  FS 

CON=l./SQRTF( 2. *3. 141592653589793) 
DO  230  1=1*1001 

230  FS ( I >=CON*EXPF( (-.5)*T(I)**2) 

-  COMPUTE  GS 


-  64  - 


n  n  n  n  n  n  n  n  c?  n  n  n  o  000000000  no  ooooooooo  no  non 


Program  E.2,  page  2 


GS(1001)=. 76198  53024  16059  16  E-23 
DO  235  1=1*91 

235  FS5( I ) =C0N*EXPF ( ( -.5 )*TS ( I ) **2 ) 

DO  250  1=2*10 
J=  11-1 
K=  10*<J-1> 

Tl  =  16067  «  * ( FSS ( K+l )  +FSS(  K+l  1 )  ) 

T2  =  106300.  * ( FSS( K+2 ) +FSS( K+10 )  ) 

T3=48525.  * ( FSS(K+3 >+FSS< K+9 )) 

T 4=272400*  *(FSS(K+4)+F$S<K+8 ) ) 

T5=260550.  * ( FSS( K+5 )+FSS( K+7 ) ) 

T 6=427368 «  *(F$S(K+6)) 

AT1=<5. *0.001) Z299376. 

AT2=T1+T2-T3+T  4-T5+T6 
A=AT1*AT2 
L  =  1 002  — I 
M=1002— I +1 
250  GS(L)=GS(M)+A 

DO  260  1=11*1001 
J=  1001-1 

Tl  =  16067.  *  ( FS  (J  +  D+FS  (J+ll)) 

T2=106300*  *  ( FS  CJ+2)+F$  (J+10)) 

T3=48525.  *<FS  (J+3J+FS  (J+9)) 

T 4=272400.  *(FS  (J+4J+FS  <J+8)J 
T5=260550.  *(FS  (J+5)+FS  (J+7H  * 

T 6=427368.  *(FS  <J+6)> 

AT1=  (-5.*0.01  )  Z299376 . 

AT2=T1+T2-T3+T4-T5+T6 
A=AT 1*AT  2 
K=1002-I 
L=1002-I+10 
260  GS(K)=GS(L)+A 

-  COMPUTE  FLS 

DO  270  1=1*1001 
270  FLS(I>*FSm-Tm*GS<  I  > 

- WRITE  FS »GS* FLS 

WRITE  OUTPUT  TAPE  6 *25* ( FS ( I ) » 1*1 *1001 ♦ 100 ) 

WRITE  OUTPUT  TAPE  6*25* ( GS< I ) *  1  =  1  *  1001 *100 ) 

WRITE  OUTPUT  TAPE  6* 25 *  (  FLSU  »  »  I  =  1 , 1001 » 100 ) 

- ROUTINE  FOR  COMPUTING  THCTT  »THCT  »FLAMCT  *FLCT  *RHOCT  *AND  RCT 

-  COMPUTE  THl  AND  ASSOCIATED  D1 

300  ASSIGN  310  TO  ISW 
DO  330  1=1*1001 

PHI  =  TU)*(GSU  )/FSU)  ) 

I F (  PHI-  .5)  330*330*305 
305  GO  TO  ISW*  (310*315) 
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310  ASSIGN  315  TO  ISW 
IFIRST  =  I 
LENGTH=1002-IFIRST 
315  TH1 ( I )  =  SQRTF  <2.*PHI-  1.) 

Dl( I )  *  THl(I)  *  T( I) 

330  CONTINUE 

WRITE  OUTPUT  TAPE  6 *10* IFIRST .LENGTH 
WRITE  OUTPUT  TAPE  6*  25, (  D1 ( I ) *  I *1 » 1001 ) 

WRITE  OUTPUT  TAPE  6*  25* < TH1 <I) »I«1 *1001 ) 

—  COMPUTE  THCTT  BYINTERPOLATION  IN  Dl«  F ( TH1 )  TABLE 

DO  450  J*I0FD1 * IOFD2  * IOFD3 
IF ( J-l ) 350*350*360 
350  THCTTIJ)*0. 

GO  TO  425 
360  DT*D(J) 

DO  390  1*1,1001 
I F(DT-D1( I ) >370*380*390 
370  TA=DT-D1 ( I-l ) 

TB=D1 ( I )-Dl( I-l) 

THCTT  ( J )*TH1 < I -1 )  +  (TA/TB)*(TH1(  D-THK  I-l)) 

GO  TO  425 

380  THCTT  < J)*TH1 ( I ) 

GO  TO  425 
390  CONTINUE 

THCTT  <  J ) =TH1 ( 1001 ) 

GO  TO  425 

425  WRITE'  OUTPUT  TAPE  6*25*THCTT  ( J  ) 

450  CONTINUE 

-  USING  THCTT  COMPUTE  THCT  BY  NEWTON-RAPHSON 

DO  1500  J*  IOFD1* IOFD2.IOFD3 
I  ND*0 

T I=D ( J) /THCTT ( J ) 

C 

1100  WRITE  OUTPUT  TAPE  6,1110*TI 
1110  FORMAT ( E20*8 ) 

L ( 1 ) =  XLOCF ( L ( 1 )  ) 

L ( 2 ) *  XLOCF ( T ( 1 )  ) 

L ( 3 ) *  XLOCF ( GS ( 1 ) ) 

L<4)=  1 
L  (  5  )  *  1 
L  ( 6 )  =  5 
L ( 7 ) *  1001 
L(8)=  0 

GOFTI  =  TABF (  TI  »L(1)  ) 

I =L I  8 ) 

GO  TO  ( 1140*1120) , I 
C 

1120  PRINT  1125*  TI 

1125  FORMAT (  6H1G  BADE20.8) 

GO  TO  9999 


-  66  - 


n  n  n 


Program  E.2,  page  4 


C 


1140  L 

(d« 

XLOCF ( 

LID) 

L 

( 2 )  * 

XLOCF ( 

Til)) 

L 

(3)» 

XLOCF! 

FS  ( 1 ) ) 

L 

( 4 )  = 

1 

L 

(5)* 

1 

L 

( 6 )  = 

5 

L 

(7)« 

1001 

L 

(8)« 

0 

FOFTI®  TABF ( TI *L(1) ) 

1=  L(8) 

GO  TO  ( 1180*1160) *  I 
C 

1160  PRINT  1165*  TI 

1165  FORMAT  (6H1F  BA0E20.8) 

GO  TO  9999 
C 

1180  QI*TI*(GOFTIVFOFTI ) 

Yl  *  QI  “  .5*! 1.  +  (D(J)/TI )**2) 
QIPR=(1.+TI**2)  *  (GOFTI/FOFTI )  -TI 
Y IPR*5  QIPR+  (  D!J)**2  /'  Tl**3) 

TINC®  TI  -  (YI/YIPR) 

C 

IF C INO-O) 1200* 1200* 1220 
1200  I ND=  1 

TI  =  TINC 
GO  TO  1100 
C 

1220  I F I  TINC-  TI 11300. 1300*1250 
C 

1250  T I =  TINC 

GO  TO  1100 
C 

1300  THCT(J)  *  D(J)/TI 

WRITE  OUTPUT  TAPE  6.1305 
1305  FORMAT { 1H  ) 

WRITE  OUTPUT  TAPE  6»1310*THCT< J) 

1310  FORMAT! E20. 8) 

C 

1500  CONTINUE 

-  COMPUTE  FLAMCT .  FLCT  USING  FS  TABLE 

DO  2500  J®  I0FD1 » I0FD2  * I0FD3 
DTH=D(J)/THCT( J) 

L ( 1 ) “XLOCF ( L ( 1 )  ) 

L  ( 2 ) ®XLOCF ( T ( 1 ) ) 

L ( 3 ) aXLOCF ( FS( 1 ) ) 

L(4)»l 
L ( 5 ) *1 
L(6)*5 
L ( 7 ) =  1001 
L(8)»  0 

FSDTH«TABF ( DTH.L ( 1 ) ) 

1=  L  ( 8 ) 

GO  TO  (2440*2420).! 

C 
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2420  PRINT 2425*J«DTH 

2425  FORMAT  (  I5H1FS ( D/THCT )  BADl4»E20.8) 

GO  TO  9999 
C 

2440  FLAMCA“ 1 • 

FLAMCB  =  ( 2.*THCT ( J )  )  /  ( (( 1 .-THCT ( J ) **2  )**2  )*  FSDTH) 
FLAMCT ( J ) *FLAMCA*FLAMCB 
IF (FLAMCT(J)-O.) 2450*2450,2445 
2445  FLCT(J)  *  L0G10F  ( FLAMCT < J I  ) 

GO  TO  2500 

2450  FLCT(J)  *  -. 99999999E+09 
C 

2500  CONTINUE 

-  COMPUTE  RHO  AND  LOGRHO 

DO  2950  J=I0FD1 * IOFD2  * IOFD3 

RHOCT ( J )*THCT ( J) **2  / ( 1 .-THCT ( J )**2  ) 

IF(RHOCT( J)-0. 12920,2920*2930 
2920  RCT ( J) *-99.9999 
GO  TO  2950 

2930  RCT ( J) =LOG10F ( RHOCT ( J  )  ) 

2950  CONTINUE 

-  PRINT  D,  THCTT, THCT, FLCT, RCT 

WRITE  OUTPUT  TAPE  6»3050»NU 
3050  FORMAT ( 4H1NU*I 4) 

WRITE  OUTPUT  TAPE  6,3100 

3100  FORMAT ( llOHO  D  THCTT 

X  THCT  FLCT  RCT) 

DO  3200  J=I0FD1 , IOFD2 , IOFD3 

WRITE  OUTPUT  TAPE  6»3150»D< J ) , THCTT ( J ) .THCT ( J ) »FLCT ( J ) ,RCT ( J ) 
3150  FORMAT (F10.1.4E25.8) 

3200  CONTINUE 

- END 

9999  CONTINUE 
CALL  EXIT 
END 
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F.  Subprogram  for  Aitkens  Interpolation  (Modification  of  SHARE 
Distribution  Number  355i408) 

FA  P 

ENTRY  TAB 
ENTRY  XTAB 

REM  XI  CONTAINS  X  WITH  ADJUSTED  SIGN 

REM  COMMON-2  IS  ERAS 

REM  K  NOT  GREATER  THAN  5 

REM  COMMON  CONTAINS  X 

REM  COMMON-1  CONTAIN  LOC  OF  LIST 

REM  COMMON  FOR  AUXILIARY 

REM  XT  ABLE  FLOATS  FROM  COMMON-8 

REM  TO  COMMON-8  WHEN  K=5 

REM  Y  TABLE  FLOATS  FROM  COMMON-9 

REM  TO  COMMON- 14 

REM  ON  OV  OR  UNDER  FLOW  OR  DCT 

REM  ORIGINAL  ARG  IN  ACC 

REM  1  LIST  (8)  IF  SOLN  IS  GOOD 

REM  2  LIST  (8)  IF  ERROR 

REM  LIST  ( 1 ) =LOCAT ION  OF  LIST  (1) 

REM  LIST  ( 2 ) =LOCATlON  OF  X  TABLE 

REM  LIST  (3)=LOCATION  OF  Y  TABLE 

REM  LIST  ( 4 ) =DEL  X 

REM  LIST  ( 5 ) =DEL  Y 

REM  LIST  ( 6 ) =K 

REM  LIST  ( 7  >  SP 

REM  LIST  (81=1-SUCCESSFUL»  2-ERROR 
J1  CLA  0*1  A*LOC  OF  X  TABLE 


J23 

J23A 

J23B 

J24 

J25 

J26 

J27 

J28 

J29 

J30 


CLA  0*1 
NOP 

CAS  XI 
TRA  J9 
TRA  J9 
TXI  J7» 1 *0 
TXI  J8»2»0 
TXL  J1»1»0 
CLA  1C  2 
LBT 

TRA  J13 
TRA  J26 
CLA  0*1 
FSB  COMMON 
SSP 

STO  COMMON-2 
T IX  J18  *1*0 
CLA  0*1 
FSB  COMMON 
SSP 

FSB  COMMON-2 
TMI  J25 
TXL  J24»l»0 
T IX  J23B.1»0 
T IX  J23 *2 *0 
TXI  J26*l *0 
T IX  J26  *2  *0 
TXI  J27 *1*0 
TXI  J28  *2*0 
LXD  XI *4 
CLA  COMMON 
FSB  0*1 


CHS  IF  XU  +  l)  LESS  THAN  XII) 

XI  HAS  SIGN  ADJUSTED 
X  LESS  THAN  X(I) 

X=X(I) 

X  GREATER  THAN  X(I)»  D*DELTA  X 
D=DELTA  Y 

KODD ( P-( K+l ) /2— 1 ) *KEVEN( P-( K/2+1 )-l ) 
LOC  OF  K  IN  ADDRESS 

K  EVEN 
K  ODD 

A* LOC  OF  XTABLE 
LOCATION  OF  X 

ERASABLE 
D3DELT a  X 
A=LOC  OF  XTAB 


X  CLOSER  TO  X<I) 

DECREMENT  SAME  AS  J8  TAB  54.1 

DECR  DELTAX  TAB  54.2 

DECR  DELTAY  TAB  54.3 

X  CLOSER  TO  X( I+1)D»DELTA  X 
D=DELTA  Y 

- ( K/2 IDELTA  X»-(K+l)/2  DELTA  X 
-(K/2JDELTA  Y.-(K+l)/2  DELTA  Y 
K1=0  *0  »  K+l 

A=LOC  OF.  X  TABLE 
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STO 

COMMON-2* 4 

J31 

TXI 

J32  ,1,0 

D=DELTA  X 

J32 

T IX 

J29»4»l 

LXD 

K 1  *4 

J34 

CLA 

0,2 

A=Y  TABLE 

STO 

COMMON- 8, 4 

J36 

TXI 

J37  ,2,0 

D*DELTA  Y 

J37 

T IX 

J34 ,4,1 

LXD 

K  » 1 

K  =  0 ,0 ,5 

J39 

DCT 

BEGIN  AITKEN  INTERPOLATION 

J40 

t6v 

J41 

J41 

PXD 

0,1 

J42 

PDX 

0,2 

J43 

TOV 

J44 

J44 

CLA 

COMMON- 2, 2 

FSB 

COMMON- 3,1 

STO 

COMMON- 2 

ERASE,  X(I)-X(I+1) 

LDQ 

COMMON- 9,1 

FMP 

COMMON- 2, 2 

STO 

K 1 

ERASE,  Y ( I )  <X-X< I+ll ) 

LDQ 

COMMON- 3,1 

FMP 

COMMON- 8, 2 

CHS 

FAD 

K 1 

TOV 

J67 

TAB 

53.1 

TOO 

J54 

TAB 

53.2 

J54 

FOP 

COMMON- 2 

TAB 

54 

TOO 

J67 

TAB 

54.1 

STQ 

COMMON-8,2 

T IX 

J43 ,2,1 

T7X 

J41»l»l 

DCT 

TRA 

J67 

DIVIDE  CHECK  ON 

CLA 

DONE 

PREPARE  NORMAL  EXIT 

J61 

STO 

0 

A=LIST  (8) 

CLA 

COMMON- 9 

LXD 

RE61, 1 

LXD 

REG2,2 

LXD 

REG4,4 

TRA 

1,4 

J67 

CLA 

COMMON 

SET  ERROR  RETURN 

STO 

COMMON- 9 

CLA 

DTWO 

TRA 

J61 

J71 

SXD 

REG 1,1 

ENTRY  POINT  FOR  TABLE 

SXD 

REG2,2 

SXD 

REG4.4 

STO 

COMMON 

X 

STQ 

COMMON- 1 

LOCATION  OF  LIST 

CLA 

COMMON- 1 

COM 

ADD 

DONE 

PDX 

0,4 

2-COMPLEMENT  OF  LIST  LOCATION 

CLA 

-1  ,4 

ADDRESS  OF  X  TABLE 

ARS 

18 

STA 

J1 

STA 

J13 
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STA  J18 
STA  J30 
STA  J137 
STA  J139 
CLA  -2*4 
ARS  IS 
STA  J34 
CLA  -5*4 
STO  K 
ADD  DONE 
STO  K1 
SUB  DONE 
ARS  18 
STA  K2 
CLA  -3*4 
STD  J6 
STD  J17 
STD  J24 
STD  J31 
STD  J23A 
CLA  -4*4 
STD  J7 
STD  J25 
STD  J36 
STD  J23B 
LXD  K»  1 
LDQ  DONE+1 » 1 
MPY  -3*4 
ALS  17 

STO  COMMON- 2 
COM  . 

ADD  DONE 
STD  J26 
LDQ  DONE+1 *1 
MPY  -4*4 
ALS  17 

STO  COMMON-8 
COM 

ADD  DONE 
STD  J27 
CLA  K2 
LBT 

TRA  J124 
TRA  J130 

J 124  CLA  COMMON- 2 
ADD  -3*4 
STO  COMMON- 2 
CLA  COMMON- 8 
ADD  -4*4 
STO  COMMON-8 
J130  CLA  -6*4 
SUB  DONE 
LRS  35 
MPY  -3*4 
ALS  17 

SUB  COMMON- 2 
STD  J8 


ACC  HAS  DELTA  X 


DELTA  Y  IN  ACC 


TAB137*1 


TAB141.1 


( AK/2 )  DELTAX 


( AK/2 JDELTAY »  TEMP  ERAS 


K  EVEN 
K  ODD 
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STD 

J23 

J 137 

CLA 

0 

LXD 

J17  » 1 

J 139 

FSB 

0.1 

TPL 

J144 

CLA 

NOP 

STO 

J2 

CLA 

COMMON 

TRA 

J 148 

J 1 44 

CLA 

CHS 

STO 

J2 

CLA 

COMMON 

CHS 

J148 

STO 

XI 

CLA 

0*4 

SUB 

SEVEN 

ARS 

18 

STA 

J61 

LXD 

COMMON-2* 1 

LXD 

COMMON- 8 *2 

CLA 

-6*4 

SUB 

K1 

TNZ 

J1 

PDX 

0*1 

PDX 

0*2 

TRA 

J28 

J156 

CLA 

LOCATE 

TRA 

1*4 

K 

PZ£ 

0*0*0 

K 1 

PZE 

0*0*0 

K2 

PZE 

0*0*0 

XI 

PZE 

PZE 

0,0.3 

PZE 

0*0*2 

DTWO 

PZE 

0*0,2 

PZE 

0*0*1 

DONE 

PZE 

0*0*1 

NOP 

NOP 

CHS 

CHS 

SEVEN 

PZE 

0*0*7 

REG1 

PZE 

REG2 

PZE 

REG4 

PZE 

LOCATE 

PZE 

0.0.J71 

TAB 

SYN 

J71 

XTAB 

SYN 

J156 

COMMON  -206 

COMMON  COMMON  9 
END 


TAB173.1 

A=XTAB, DETERMINE  DIRECTION 
OF  MONOTONICITY 
X ( 1 ) -X ( 2  ) 


TAB154.1 
TAB154.2 
TAB154.3 
TAB154.4 
TAB154.5 
TAB  155 


LOCATING  SUB  ROUTINE 
0  *0  ,  K 
0*0*  K+l 
K  *0  *  0 

X  WITH  ADJUSTED  SIGN 
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G. 1  Optimal  Sample  Size  -  First  Search. Procedure 


-  OPTIMAL  5AMPLE  SIZE  -  FIRST  SEARCH  PROCEDURE 

5  FORMAT (110) 

10  FORMAT (F5*1*E30*20»2F5«2»I5  ) 

20  FORMAT ( E40.20 ) 

30  FORMAT ( 1H  ) 

C 

READ  5 »NUMB 
C 

DO  800  IMAIN=1 »NUMB 
C 

100  READ  10  »FLNU*BETA»D»FL»MAXT  T 
C 

A=-4»*( FLNU+D*#2 ) 

B=-U6.-FLNU)*D**2+FLNU) 

C=D**2* (FLNU-2 • ) 

RHOLO* ( -B-SQRT F ( B**2- ( 4 . *A*C >)) / ( 2 • *A ) 
FLAM=10.**FL 

CON=(FLNU**< «5*FLNU) )/8ETA 
EXP=-.5*(FLNU+1. ) 

FSD=CON*(FLNU+D**2)**EXP 
RHOHl*SQRTF( .5*FLAM*FSD) 

PRINT  20»A*B*C» RHOLO ♦ FL AM » CON » EXP » FSD » RHOH I 

C 

DO  500  I Ts 1» MAX  IT 
I F ( IT— 2 )110»120»130 
110  RHO*RHOLO 
GO  TO  150 
120  RHO=RHOHI 
GO  TO  150 
130  RHO=RHOT 
GO  TO  150 
C 

150  THETA=SQRTF(RHO/(RHO+l.» ) 

FSDTH*CON*(FLNU+< D/THETA )**2)**EXP 
Tl= ( FLNU+T D/THETA) **2 ) *FSDTH/ ( FLNU-1 . ) 

T2=l »/SQRTF ( RHO* ( RHO+1 • ) **3 ) 

GPR= ( .5*FLAM*T1*T2 )-l . 

PRINT  20»THETA»FSDTH  *T1»T2*GPR 
C 

I F ( I T-2 ) 210  *  220 *300 
210  RHOL*RHOLO 
GPRL=GPR 
GO  TO  500 
220  RHOH=RHOHl 
GPRH=GPR 
GO  TO  400 
C 

300  IFIGPR-O* ) 310*320*330 
310  RHOH=RHOT 
GPRH=GPR 
GO  TO  400 
320  RHOOPaRHOT 

PRINT  20 *RHOOP *GPR 
GO  TO  800 
330  RHOL=RHOT 
GPRL=GPR 
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GO  TO  400 
C 

400  TOP* ( GPRL*RHOH ) - ( GPRH*RHOL ) 
BOT=GPRL-GPRH 
RHOT*TOP/BOT 
PRINT  20*TOP*BOT»RHOT 

L  RL=« 43429448 190325 182765*LOGF ( RHOL ) 

RH«.4342944B190325182765*LOGF(RHOH) 
RT*. 43429448 190325 182765*LOGF  <  RHOT  ) 
C 

PRINT  20 »RHOL »GPRL »RHOH,GPRH»RHOT 
PRINT  20*RL»RH»RT 
PRINT  30 
C 

500  CONTINUE 
C 

800  CONTINUE 
C 

END 
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G .2  Optimal  Sample  Size  -  Second  Search  Procedure 


-  OPTIMAL  SAMPLE  SIZE  -  SECOND  SEARCH  PROCEDURE 

5  FORMAT (110) 

10  FORMAT (F5*1*E30«20*2F5«2»I5) 

20  FORMAT (3E40. 20) 

30  FORMAT ( 1H  ) 

C 

READ  5  * NUMB 
C 

DO  1000  IMAIN=1 »NUMB 
C 

READ  10*FLNU»BETA»D»FL*MAXIT 
PRINT10*FLNU»BETA»D»FL»MAXI T 
C 

A=-4.*(FLNU+D**2) 

B*-( (6.-FLNU)*D**2+FLNU) 

OD**2*(FLNU-2.) 

RHOLO»(-B-SQRTF<B**2-(4.*A*C) >  >/(2.*A) 
FLAM=10.**FL 

CON= ( FLNU** ( • 5*FLNU ) ) /BETA 
EXP=-,5*(FLNU+1. ) 

FSD=CON* ( FLNU+D**2 ) **EXP 

RHOH I =SQRTF ( . 5*FLAM* ( FLNU+D**2 ) *FSD/ ( FLNU-1 . ) ) 
PRINT  20  »A*B»C  »RHOLO»FLAM*CON»EXP*FSD»RHOHI 

C 

1310=0 

1330*0 

DO  500  IT=1»MAXIT 
IF ( IT— 2 )110» 12 0*130 
110  RHO=RHOLO 
GO  TO  150 
120  RHO=RHOHI 
GO  TO  150 
130  RHO=RHOT 
GO  TO  150 
C 

150  THETA*SQRTF(RHO/(RHO+l.) ) 

FSDTH»CON* ( FLNU+ ( D/THETA ) **2 ) **EXP 
T 1* ( FLNU+ ( D/THETA ) **2 )*FSDTH/ ( FLNU-1 • ) 
T2=l./SQRTF(RHO*(RHO+l. )**3) 

GPR= ( «  5*FLAM*T1*T2 ) -1 • 

PRINT  20*GPR 
C 

IF( IT-2)210»220»300 
210  RHOL=RHOLO 

I F (GPR— 0* J212#212*214 
212  RHOT*0. 

GO  TO  900 
214  GPRL=GPR 
GPLB4=GPR 
GO  TO  500 
220  RHOH=RHOHl 
GPRH*6PR 
6PHB4=GPR 
GO  TO  400 
C 

300  I F ( GPR-0* ) 310*400*330 
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310  RH0H*RH0T 
GPRH-GPR 

IF (GPRH-GPHB4) 900.900.315 
315  GPHB4=GPRH 
1310=1310+1 
I 330*0 

I F( 1310-3)400*450.450 
330  RHOL=RHOT 
GPRL=GPR 

I F(GPLB4-GPRL) 900*900. 335 
335  GPLB4*GPRL 
1330*1330+1 
1310*0 

I F ( 1330-3)400*450.450 
C 

400  RH0T=( ( GPRL*RH0H ) - ( GPRH*RHOL )  )  /  ( GPRL-GPRH ) 
PRINT  20»RHOL*RHOH*RHOT 
PRINT  30 
GO  TO  900 
C 

450  RHOT  *• 5*( RHOH+RHOL ) 

PRINT  20.RHOL.RHOH.RHOT 
PRINT  30 
1310=0 
1330=0 
GO  TO  500 
C 

500  CONTINUE 

C  ' 

900  RT **43429448 190325 18276 5*LOGF (RHOT  V 
PRINT  20 *GPR 
PRINT  20*RHOL»RHOH 
PRINT  20 *RHOT »RT 
PRINT  5  » IT 
C 

1000  CONTINUE 
END 
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-  OPTIMAL  SAMPLE  SIZE  -  NEWTON-RAPHSON  PROCEDURE 

10  FORMAT (2E30.20) 

15  FORMAT < 1H  ) 

C 

READ  10»FLNU.BETA 
PRINT10»FLNU».BETA 
READ  10  »D»FL 
PRINT10*D»FL 
READ  10  *RCT 
PRINTlOtRCT 
RHONEW=10.**RCT 
PRINT  10.RHONEW 
C 

DO  500  IMAIN=1 ♦ 15 
RHO=RHONEW 
C 

FLAM=10 •**FL 

THETA=SQRTF(RHO/(l.+RHO) ) 

CON= 1» / ( SQRTF ( FLNU ) *BET  A ) 

EXP=(-«5)*(  FLNU+1 • ) 

DTH*D/THET  A 

FSDTH=CON  *  (l.+DTH**2/FLNU)**EXP 
TA=(.5*FLAM)  /  (SQRTF (RHO*(RHO+l* 

TB* ( FLNU+DTH**2 )  /  (FLNU-1.) 

6PR= ( TA*T8*FSDTH ) -1 • 

TA=FLAM/  (4.*(FLNU-1.) ) 

TB*1»/SQRTF( (RHO*(RHO+l. >  >**5. ) 

TQA=-4. * ( FLNU+D**2 ) *RHO**2 
TQB=(D**2*FLNU-6.*D**2-FLNU)  *  RHO 
TQC“D**2*( FLNU— 2 • ) 

G2PR=TA*TB*FSDTH*  ( TQA+TUB+TQC ) 

PRINT  10»FLAM*THETA*CON*EXP  *DTH»FSDTH»TA*TB*TQA  *TQB*TOC»G2PR 
C 

RHONEW=RHO  -  (GPR/G2PR) 

RCTNEW=. 43429  44819  03251  82765  *  LOGF ( RHONEW ) 

C 

PRINT  10 »GPR  *G2PR 
PRINT  10*RH0NEWtRCTNEW 
PRINT  15 
C 

500  CONTINUE 
END 
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G.h  Net  Gain  for  Given  Optimal  Sample  Size.  Real  V  <  00  .  t  <  °° 

C  -  NET  GAIN  FOR  GIVEN  OPTIMAL  SAMPLE  SIZE*  REAL  NU.  CONVERGES  AS 

C  -  POWER  SERIES  IN  NU/(NU+T**2> 

C 

10  FORMAT (1 10) 

17  FORMAT (2E30»2Q*F10*4) 

18  F0RMAT(2E30. 20*15) 

19  FORMAT ( 1H  ) 

C 

READ  17»TNU»TBETA 
PRINT17*TNU»TBETA 
READ  17*TD*TL 
PRINT17*TD*TL 
PRINT  19 

1  READ  17*  TRHO 
PRINT  17  *TRHO 

C 

C - TFLAM.TTHET  *  T 

C 

TFLAM*10***TL 

TTHET*SQRTF( TRHO/ ( TRHO+1* ) ) 

T  *TD/TTHET 
C 

C  -  DENSITY  FS 

C 

Z*(1./SQRTF(TNU) )*(1./TBETA) 

EXP«(-.5)*(TNU+1.) 

FS*Z*<1.+T**2/TNU)**EXP 
PRINT  17  *FS 
C 

C  -  RIGHT  TAIL  GS 

C 

EPS=.1E-10  *  10.**<-T) 

C 

FRAC*T**2/TNU 
CON*  SQRTF(FRAC)  /  TBETA 
DELTA*  EPS/CON 
SUM=1 
TERM=1 
C 

DO  240  J*2  *1000 
F  J=J 

TINCR*  -U2.*FJ-3,  )/(2.*FJ-l.» )  *  (l./(FJ-l.)) 

X  *(.5*<TNU+2.*FJ-3.)» 

TERM*  TERM*T I NCR*FR AC 
SUM*TERM+SUM 
ERFAC*ABSF (TERM) 

IF (ERFAC  -DELTA)  245*240*240 
240  CONTINUE 
GO  TO  999 
C 

245  GS=.5-C0N*SUM 
C 

PRINT  18»GS,EPS»J 
C 

C  -  LINEAR  LOSS  FLS 

C 

FLS* ( ( TNU+T»*2 ) / ( TNU~1» ) ) *FS  -  T*GS 
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PRINT  17  *FLS 
PRINT  19 

-  NET  GAIN 

FNG* ( TFLAM*TTHET*FLS)-TRHO 
PRINT  17  »FNG 

FG*. 43429  44819  03251  82765  *  LOGF(FNG) 
PRINT  17.FG 
PRINT  19 

999  CONTINUE 
GO  TO  1 
END 
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G.5  Net  Gain  for  Given  Optimal  Sample  Size.  Real  v  <  «  .  t  <  v 

-  NET  GAIN  FOR  GIVEN  OPTIMAL  SAMPLE  SIZE»  REAL  NU.  CONVERGES  FOR 

-  T**2  LESS  THAN  NU 

10  FORMAT (110) 

17  FORMAT (2E30.20.F10.3) 

18  FORMATI2E30. 20*15) 

19  FORMAT ( 1H  ) 

READ  17*TNU*TBETA 
PRINT17»TNU*TBETA 
REAO  17  *TD»TL 
PRINT17.TD.TL 
PRINT  19 

1  READ  17»TRHO 
PR INT17  *TRHO 

-  TFLAM.TTHET.T 

TFLAM*10.**TL 

TTHET-SQRTF  < TRHO/ ( TRHO+1 •) ) 

T  *TD/TTHET 

-  DENSITY  FS 

Z*(1./SQRTF(TNU))*(1./TBETA> 

EXP*(-.5)*(TNU+1.) 

FS*Z*( l.+T**2/TNU)**EXP 
PRINT  17  *FS 

-  RIGHT  TAIL  GS 

EPS=.1E— 10  *  10.**<-T> 

FRAC*  TNU/(TNU+T**2) 

CON  *  (FRAC**(.5*TNU) I/TBETA 
DELTA  *  <T**2/TNU>  *  (EPS/CON) 

ERFAC*1. 

SUM*  l./TNU 
TERM*  l./TNU 

DO  340  J*2 ♦ 1000 
FJ*J 

T INCR*  .5*(  C TNU+2.*FJ-4. ) / ( TNU+2.*FJ-2. ) ) 

X  *(<2.*FJ-3.)/(FJ-l.n 

TERM*  TERM*T INCR*FRAC 
SUM*  SUM+TERM 
ERFAC*  ERFAC*FRAC 
IF(ERFAC-DELTA)  345.340,340 
340  CONTINUE 
GO  TO  999 

345  GS*C0N*SUM 

PRINT  18.GS.EPS.J 

-  LINEAR  LOSS  FLS 


-  8o  - 


non 


Program  G.5>  page  2 


FLS* I ( TNU+T**2 ) / ( TNU-1 • ) )*FS  -  T*GS 
PRINT  17.FLS 
PRINT  19 

-  NET  GAIN 

FNG«(TFLAM*TTHET#FLS)-TRHO 
FG=. 43429  44819  03251  82765  *  LOGF ( FNGJ 
PRINT  17»FNG.FG 
PRINT  19 

999  CONTINUE 
GO  TO  1 
END 
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